My Project
PreconditionerFactory.hpp
1 
2 /*
3  Copyright 2009, 2010 SINTEF ICT, Applied Mathematics.
4  Copyright 2019 SINTEF Digital, Mathematics and Cybernetics.
5 
6  This file is part of the Open Porous Media project (OPM).
7 
8  OPM is free software: you can redistribute it and/or modify
9  it under the terms of the GNU General Public License as published by
10  the Free Software Foundation, either version 3 of the License, or
11  (at your option) any later version.
12 
13  OPM is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU General Public License for more details.
17 
18  You should have received a copy of the GNU General Public License
19  along with OPM. If not, see <http://www.gnu.org/licenses/>.
20 */
21 
22 #ifndef OPM_PRECONDITIONERFACTORY_HEADER
23 #define OPM_PRECONDITIONERFACTORY_HEADER
24 
25 #include <opm/simulators/linalg/OwningBlockPreconditioner.hpp>
26 #include <opm/simulators/linalg/OwningTwoLevelPreconditioner.hpp>
27 #include <opm/simulators/linalg/ParallelOverlappingILU0.hpp>
28 #include <opm/simulators/linalg/PreconditionerWithUpdate.hpp>
29 #include <opm/simulators/linalg/PropertyTree.hpp>
31 
32 #include <dune/istl/paamg/amg.hh>
33 #include <dune/istl/paamg/kamg.hh>
34 #include <dune/istl/paamg/fastamg.hh>
35 #include <dune/istl/preconditioners.hh>
36 
37 #include <map>
38 #include <memory>
39 #include <limits>
40 
41 namespace Opm
42 {
43 
49 template <class Operator, class Comm>
51 {
52 public:
54  using Matrix = typename Operator::matrix_type;
55  using Vector = typename Operator::domain_type; // Assuming symmetry: that domain and range types are the same.
56 
58  using PrecPtr = std::shared_ptr<Dune::PreconditionerWithUpdate<Vector, Vector>>;
59 
61  using Creator = std::function<PrecPtr(const Operator&, const PropertyTree&,
62  const std::function<Vector()>&, std::size_t)>;
63  using ParCreator = std::function<PrecPtr(const Operator&, const PropertyTree&,
64  const std::function<Vector()>&, std::size_t, const Comm&)>;
65 
71  static PrecPtr create(const Operator& op, const PropertyTree& prm,
72  const std::function<Vector()>& weightsCalculator = {},
73  std::size_t pressureIndex = std::numeric_limits<std::size_t>::max())
74  {
75  return instance().doCreate(op, prm, weightsCalculator, pressureIndex);
76  }
77 
84  static PrecPtr create(const Operator& op, const PropertyTree& prm,
85  const std::function<Vector()>& weightsCalculator, const Comm& comm,
86  std::size_t pressureIndex = std::numeric_limits<std::size_t>::max())
87  {
88  return instance().doCreate(op, prm, weightsCalculator, pressureIndex, comm);
89  }
90 
96  static PrecPtr create(const Operator& op, const PropertyTree& prm, const Comm& comm,
97  std::size_t pressureIndex = std::numeric_limits<std::size_t>::max())
98  {
99  return instance().doCreate(op, prm, std::function<Vector()>(), pressureIndex, comm);
100  }
108  static void addCreator(const std::string& type, Creator creator)
109  {
110  instance().doAddCreator(type, creator);
111  }
112 
120  static void addCreator(const std::string& type, ParCreator creator)
121  {
122  instance().doAddCreator(type, creator);
123  }
124 
125  using CriterionBase
126  = Dune::Amg::AggregationCriterion<Dune::Amg::SymmetricDependency<Matrix, Dune::Amg::FirstDiagonal>>;
127  using Criterion = Dune::Amg::CoarsenCriterion<CriterionBase>;
128 
129  // Helpers for creation of AMG preconditioner.
130  static Criterion amgCriterion(const PropertyTree& prm)
131  {
132  Criterion criterion(15, prm.get<int>("coarsenTarget", 1200));
133  criterion.setDefaultValuesIsotropic(2);
134  criterion.setAlpha(prm.get<double>("alpha", 0.33));
135  criterion.setBeta(prm.get<double>("beta", 1e-5));
136  criterion.setMaxLevel(prm.get<int>("maxlevel", 15));
137  criterion.setSkipIsolated(prm.get<bool>("skip_isolated", false));
138  criterion.setNoPreSmoothSteps(prm.get<int>("pre_smooth", 1));
139  criterion.setNoPostSmoothSteps(prm.get<int>("post_smooth", 1));
140  criterion.setDebugLevel(prm.get<int>("verbosity", 0));
141  // As the default we request to accumulate data to 1 process always as our matrix
142  // graph might be unsymmetric and hence not supported by the PTScotch/ParMetis
143  // calls in DUNE. Accumulating to 1 skips PTScotch/ParMetis
144  criterion.setAccumulate(static_cast<Dune::Amg::AccumulationMode>(prm.get<int>("accumulate", 1)));
145  criterion.setProlongationDampingFactor(prm.get<double>("prolongationdamping", 1.6));
146  criterion.setMaxDistance(prm.get<int>("maxdistance", 2));
147  criterion.setMaxConnectivity(prm.get<int>("maxconnectivity", 15));
148  criterion.setMaxAggregateSize(prm.get<int>("maxaggsize", 6));
149  criterion.setMinAggregateSize(prm.get<int>("minaggsize", 4));
150  return criterion;
151  }
152 private:
153 
156  template <typename X> struct Id { using Type = X; };
157 
158  template <typename Smoother>
159  static auto amgSmootherArgs(const PropertyTree& prm)
160  {
161  return amgSmootherArgs(prm, Id<Smoother>());
162  }
163 
164  template <typename Smoother>
165  static auto amgSmootherArgs(const PropertyTree& prm,
166  Id<Smoother>)
167  {
168  using SmootherArgs = typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
169  SmootherArgs smootherArgs;
170  smootherArgs.iterations = prm.get<int>("iterations", 1);
171  // smootherArgs.overlap=SmootherArgs::vertex;
172  // smootherArgs.overlap=SmootherArgs::none;
173  // smootherArgs.overlap=SmootherArgs::aggregate;
174  smootherArgs.relaxationFactor = prm.get<double>("relaxation", 1.0);
175  return smootherArgs;
176  }
177 
178  static auto amgSmootherArgs(const PropertyTree& prm,
180  {
182  using SmootherArgs = typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
183  SmootherArgs smootherArgs;
184  smootherArgs.iterations = prm.get<int>("iterations", 1);
185  const int iluwitdh = prm.get<int>("iluwidth", 0);
186  smootherArgs.setN(iluwitdh);
187  const MILU_VARIANT milu = convertString2Milu(prm.get<std::string>("milutype", std::string("ilu")));
188  smootherArgs.setMilu(milu);
189  // smootherArgs.overlap=SmootherArgs::vertex;
190  // smootherArgs.overlap=SmootherArgs::none;
191  // smootherArgs.overlap=SmootherArgs::aggregate;
192  smootherArgs.relaxationFactor = prm.get<double>("relaxation", 1.0);
193  return smootherArgs;
194  }
195 
196  template <class Smoother>
197  static PrecPtr makeAmgPreconditioner(const Operator& op, const PropertyTree& prm, bool useKamg = false)
198  {
199  auto crit = amgCriterion(prm);
200  auto sargs = amgSmootherArgs<Smoother>(prm);
201  if(useKamg){
202  return std::make_shared<
205  >
206  >(op, crit, sargs,
207  prm.get<size_t>("max_krylov", 1),
208  prm.get<double>("min_reduction", 1e-1) );
209  }else{
210  return std::make_shared<Dune::Amg::AMGCPR<Operator, Vector, Smoother>>(op, crit, sargs);
211  }
212  }
213 
218  template <class CommArg>
219  static size_t interiorIfGhostLast(const CommArg& comm)
220  {
221  size_t interior_count = 0;
222  size_t highest_interior_index = 0;
223  const auto& is = comm.indexSet();
224  for (const auto& ind : is) {
225  if (CommArg::OwnerSet::contains(ind.local().attribute())) {
226  ++interior_count;
227  highest_interior_index = std::max(highest_interior_index, ind.local().local());
228  }
229  }
230  if (highest_interior_index + 1 == interior_count) {
231  return interior_count;
232  } else {
233  return is.size();
234  }
235  }
236 
237  static PrecPtr
238  createParILU(const Operator& op, const PropertyTree& prm, const Comm& comm, const int ilulevel)
239  {
240  const double w = prm.get<double>("relaxation", 1.0);
241  const bool redblack = prm.get<bool>("redblack", false);
242  const bool reorder_spheres = prm.get<bool>("reorder_spheres", false);
243  // Already a parallel preconditioner. Need to pass comm, but no need to wrap it in a BlockPreconditioner.
244  if (ilulevel == 0) {
245  const size_t num_interior = interiorIfGhostLast(comm);
246  return std::make_shared<Opm::ParallelOverlappingILU0<Matrix, Vector, Vector, Comm>>(
247  op.getmat(), comm, w, Opm::MILU_VARIANT::ILU, num_interior, redblack, reorder_spheres);
248  } else {
249  return std::make_shared<Opm::ParallelOverlappingILU0<Matrix, Vector, Vector, Comm>>(
250  op.getmat(), comm, ilulevel, w, Opm::MILU_VARIANT::ILU, redblack, reorder_spheres);
251  }
252  }
253 
254  // Add a useful default set of preconditioners to the factory.
255  // This is the default template, used for parallel preconditioners.
256  // (Serial specialization below).
257  template <class CommArg>
258  void addStandardPreconditioners(const CommArg*)
259  {
260  using namespace Dune;
261  using O = Operator;
262  using M = Matrix;
263  using V = Vector;
264  using P = PropertyTree;
265  using C = Comm;
266  doAddCreator("ILU0", [](const O& op, const P& prm, const std::function<Vector()>&, std::size_t, const C& comm) {
267  return createParILU(op, prm, comm, 0);
268  });
269  doAddCreator("ParOverILU0", [](const O& op, const P& prm, const std::function<Vector()>&, std::size_t, const C& comm) {
270  return createParILU(op, prm, comm, prm.get<int>("ilulevel", 0));
271  });
272  doAddCreator("ILUn", [](const O& op, const P& prm, const std::function<Vector()>&, std::size_t, const C& comm) {
273  return createParILU(op, prm, comm, prm.get<int>("ilulevel", 0));
274  });
275  doAddCreator("Jac", [](const O& op, const P& prm, const std::function<Vector()>&,
276  std::size_t, const C& comm) {
277  const int n = prm.get<int>("repeats", 1);
278  const double w = prm.get<double>("relaxation", 1.0);
279  return wrapBlockPreconditioner<DummyUpdatePreconditioner<SeqJac<M, V, V>>>(comm, op.getmat(), n, w);
280  });
281  doAddCreator("GS", [](const O& op, const P& prm, const std::function<Vector()>&, std::size_t, const C& comm) {
282  const int n = prm.get<int>("repeats", 1);
283  const double w = prm.get<double>("relaxation", 1.0);
284  return wrapBlockPreconditioner<DummyUpdatePreconditioner<SeqGS<M, V, V>>>(comm, op.getmat(), n, w);
285  });
286  doAddCreator("SOR", [](const O& op, const P& prm, const std::function<Vector()>&, std::size_t, const C& comm) {
287  const int n = prm.get<int>("repeats", 1);
288  const double w = prm.get<double>("relaxation", 1.0);
289  return wrapBlockPreconditioner<DummyUpdatePreconditioner<SeqSOR<M, V, V>>>(comm, op.getmat(), n, w);
290  });
291  doAddCreator("SSOR", [](const O& op, const P& prm, const std::function<Vector()>&, std::size_t, const C& comm) {
292  const int n = prm.get<int>("repeats", 1);
293  const double w = prm.get<double>("relaxation", 1.0);
294  return wrapBlockPreconditioner<DummyUpdatePreconditioner<SeqSSOR<M, V, V>>>(comm, op.getmat(), n, w);
295  });
296 
297  // Only add AMG preconditioners to the factory if the operator
298  // is the overlapping schwarz operator. This could be extended
299  // later, but at this point no other operators are compatible
300  // with the AMG hierarchy construction.
301  if constexpr (std::is_same_v<O, Dune::OverlappingSchwarzOperator<M, V, V, C>>) {
302  doAddCreator("amg", [](const O& op, const P& prm, const std::function<Vector()>&, std::size_t, const C& comm) {
303  const std::string smoother = prm.get<std::string>("smoother", "ParOverILU0");
304  if (smoother == "ILU0" || smoother == "ParOverILU0") {
306  auto crit = amgCriterion(prm);
307  auto sargs = amgSmootherArgs<Smoother>(prm);
308  return std::make_shared<Dune::Amg::AMGCPR<O, V, Smoother, C>>(op, crit, sargs, comm);
309  } else {
310  OPM_THROW(std::invalid_argument, "Properties: No smoother with name " << smoother << ".");
311  }
312  });
313  }
314 
315  doAddCreator("cpr", [](const O& op, const P& prm, const std::function<Vector()> weightsCalculator, std::size_t pressureIndex, const C& comm) {
316  assert(weightsCalculator);
317  if (pressureIndex == std::numeric_limits<std::size_t>::max())
318  {
319  OPM_THROW(std::logic_error, "Pressure index out of bounds. It needs to specified for CPR");
320  }
321  return std::make_shared<OwningTwoLevelPreconditioner<O, V, false, Comm>>(op, prm, weightsCalculator, pressureIndex, comm);
322  });
323  doAddCreator("cprt", [](const O& op, const P& prm, const std::function<Vector()> weightsCalculator, std::size_t pressureIndex, const C& comm) {
324  assert(weightsCalculator);
325  if (pressureIndex == std::numeric_limits<std::size_t>::max())
326  {
327  OPM_THROW(std::logic_error, "Pressure index out of bounds. It needs to specified for CPR");
328  }
329  return std::make_shared<OwningTwoLevelPreconditioner<O, V, true, Comm>>(op, prm, weightsCalculator, pressureIndex, comm);
330  });
331  }
332 
333  // Add a useful default set of preconditioners to the factory.
334  // This is the specialization for the serial case.
335  void addStandardPreconditioners(const Dune::Amg::SequentialInformation*)
336  {
337  using namespace Dune;
338  using O = Operator;
339  using M = Matrix;
340  using V = Vector;
341  using P = PropertyTree;
342  doAddCreator("ILU0", [](const O& op, const P& prm, const std::function<Vector()>&, std::size_t) {
343  const double w = prm.get<double>("relaxation", 1.0);
344  return std::make_shared<Opm::ParallelOverlappingILU0<M, V, V>>(
345  op.getmat(), 0, w, Opm::MILU_VARIANT::ILU);
346  });
347  doAddCreator("ParOverILU0", [](const O& op, const P& prm, const std::function<Vector()>&, std::size_t) {
348  const double w = prm.get<double>("relaxation", 1.0);
349  const int n = prm.get<int>("ilulevel", 0);
350  return std::make_shared<Opm::ParallelOverlappingILU0<M, V, V>>(
351  op.getmat(), n, w, Opm::MILU_VARIANT::ILU);
352  });
353  doAddCreator("ILUn", [](const O& op, const P& prm, const std::function<Vector()>&, std::size_t) {
354  const int n = prm.get<int>("ilulevel", 0);
355  const double w = prm.get<double>("relaxation", 1.0);
356  return std::make_shared<Opm::ParallelOverlappingILU0<M, V, V>>(
357  op.getmat(), n, w, Opm::MILU_VARIANT::ILU);
358  });
359  doAddCreator("Jac", [](const O& op, const P& prm, const std::function<Vector()>&, std::size_t) {
360  const int n = prm.get<int>("repeats", 1);
361  const double w = prm.get<double>("relaxation", 1.0);
362  return wrapPreconditioner<SeqJac<M, V, V>>(op.getmat(), n, w);
363  });
364  doAddCreator("GS", [](const O& op, const P& prm, const std::function<Vector()>&, std::size_t) {
365  const int n = prm.get<int>("repeats", 1);
366  const double w = prm.get<double>("relaxation", 1.0);
367  return wrapPreconditioner<SeqGS<M, V, V>>(op.getmat(), n, w);
368  });
369  doAddCreator("SOR", [](const O& op, const P& prm, const std::function<Vector()>&, std::size_t) {
370  const int n = prm.get<int>("repeats", 1);
371  const double w = prm.get<double>("relaxation", 1.0);
372  return wrapPreconditioner<SeqSOR<M, V, V>>(op.getmat(), n, w);
373  });
374  doAddCreator("SSOR", [](const O& op, const P& prm, const std::function<Vector()>&, std::size_t) {
375  const int n = prm.get<int>("repeats", 1);
376  const double w = prm.get<double>("relaxation", 1.0);
377  return wrapPreconditioner<SeqSSOR<M, V, V>>(op.getmat(), n, w);
378  });
379 
380  // Only add AMG preconditioners to the factory if the operator
381  // is an actual matrix operator.
382  if constexpr (std::is_same_v<O, Dune::MatrixAdapter<M, V, V>>) {
383  doAddCreator("amg", [](const O& op, const P& prm, const std::function<Vector()>&, std::size_t) {
384  const std::string smoother = prm.get<std::string>("smoother", "ParOverILU0");
385  if (smoother == "ILU0" || smoother == "ParOverILU0") {
386 #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
387  using Smoother = SeqILU<M, V, V>;
388 #else
389  using Smoother = SeqILU0<M, V, V>;
390 #endif
391  return makeAmgPreconditioner<Smoother>(op, prm);
392  } else if (smoother == "Jac") {
393  using Smoother = SeqJac<M, V, V>;
394  return makeAmgPreconditioner<Smoother>(op, prm);
395  } else if (smoother == "SOR") {
396  using Smoother = SeqSOR<M, V, V>;
397  return makeAmgPreconditioner<Smoother>(op, prm);
398  } else if (smoother == "SSOR") {
399  using Smoother = SeqSSOR<M, V, V>;
400  return makeAmgPreconditioner<Smoother>(op, prm);
401  } else if (smoother == "ILUn") {
402 #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
403  using Smoother = SeqILU<M, V, V>;
404 #else
405  using Smoother = SeqILUn<M, V, V>;
406 #endif
407  return makeAmgPreconditioner<Smoother>(op, prm);
408  } else {
409  OPM_THROW(std::invalid_argument, "Properties: No smoother with name " << smoother << ".");
410  }
411  });
412  doAddCreator("kamg", [](const O& op, const P& prm, const std::function<Vector()>&, std::size_t) {
413  const std::string smoother = prm.get<std::string>("smoother", "ParOverILU0");
414  if (smoother == "ILU0" || smoother == "ParOverILU0") {
415 #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
416  using Smoother = SeqILU<M, V, V>;
417 #else
418  using Smoother = SeqILU0<M, V, V>;
419 #endif
420  return makeAmgPreconditioner<Smoother>(op, prm, true);
421  } else if (smoother == "Jac") {
422  using Smoother = SeqJac<M, V, V>;
423  return makeAmgPreconditioner<Smoother>(op, prm, true);
424  } else if (smoother == "SOR") {
425  using Smoother = SeqSOR<M, V, V>;
426  return makeAmgPreconditioner<Smoother>(op, prm, true);
427  // } else if (smoother == "GS") {
428  // using Smoother = SeqGS<M, V, V>;
429  // return makeAmgPreconditioner<Smoother>(op, prm, true);
430  } else if (smoother == "SSOR") {
431  using Smoother = SeqSSOR<M, V, V>;
432  return makeAmgPreconditioner<Smoother>(op, prm, true);
433  } else if (smoother == "ILUn") {
434 #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
435  using Smoother = SeqILU<M, V, V>;
436 #else
437  using Smoother = SeqILUn<M, V, V>;
438 #endif
439  return makeAmgPreconditioner<Smoother>(op, prm, true);
440  } else {
441  OPM_THROW(std::invalid_argument, "Properties: No smoother with name " << smoother << ".");
442  }
443  });
444  doAddCreator("famg", [](const O& op, const P& prm, const std::function<Vector()>&, std::size_t) {
445  auto crit = amgCriterion(prm);
446  Dune::Amg::Parameters parms;
447  parms.setNoPreSmoothSteps(1);
448  parms.setNoPostSmoothSteps(1);
449  return wrapPreconditioner<Dune::Amg::FastAMG<O, V>>(op, crit, parms);
450  });
451  }
452  doAddCreator("cpr", [](const O& op, const P& prm, const std::function<Vector()>& weightsCalculator, std::size_t pressureIndex) {
453  if (pressureIndex == std::numeric_limits<std::size_t>::max())
454  {
455  OPM_THROW(std::logic_error, "Pressure index out of bounds. It needs to specified for CPR");
456  }
457  return std::make_shared<OwningTwoLevelPreconditioner<O, V, false>>(op, prm, weightsCalculator, pressureIndex);
458  });
459  doAddCreator("cprt", [](const O& op, const P& prm, const std::function<Vector()>& weightsCalculator, std::size_t pressureIndex) {
460  if (pressureIndex == std::numeric_limits<std::size_t>::max())
461  {
462  OPM_THROW(std::logic_error, "Pressure index out of bounds. It needs to specified for CPR");
463  }
464  return std::make_shared<OwningTwoLevelPreconditioner<O, V, true>>(op, prm, weightsCalculator, pressureIndex);
465  });
466  }
467 
468 
469  // The method that implements the singleton pattern,
470  // using the Meyers singleton technique.
471  static PreconditionerFactory& instance()
472  {
473  static PreconditionerFactory singleton;
474  return singleton;
475  }
476 
477  // Private constructor, to keep users from creating a PreconditionerFactory.
478  PreconditionerFactory()
479  {
480  Comm* dummy = nullptr;
481  addStandardPreconditioners(dummy);
482  }
483 
484  // Actually creates the product object.
485  PrecPtr doCreate(const Operator& op, const PropertyTree& prm,
486  const std::function<Vector()> weightsCalculator,
487  std::size_t pressureIndex)
488  {
489  const std::string& type = prm.get<std::string>("type", "ParOverILU0");
490  auto it = creators_.find(type);
491  if (it == creators_.end()) {
492  std::ostringstream msg;
493  msg << "Preconditioner type " << type << " is not registered in the factory. Available types are: ";
494  for (const auto& prec : creators_) {
495  msg << prec.first << ' ';
496  }
497  msg << std::endl;
498  OPM_THROW(std::invalid_argument, msg.str());
499  }
500  return it->second(op, prm, weightsCalculator, pressureIndex);
501  }
502 
503  PrecPtr doCreate(const Operator& op, const PropertyTree& prm,
504  const std::function<Vector()> weightsCalculator,
505  std::size_t pressureIndex, const Comm& comm)
506  {
507  const std::string& type = prm.get<std::string>("type", "ParOverILU0");
508  auto it = parallel_creators_.find(type);
509  if (it == parallel_creators_.end()) {
510  std::ostringstream msg;
511  msg << "Parallel preconditioner type " << type
512  << " is not registered in the factory. Available types are: ";
513  for (const auto& prec : parallel_creators_) {
514  msg << prec.first << ' ';
515  }
516  msg << std::endl;
517  OPM_THROW(std::invalid_argument, msg.str());
518  }
519  return it->second(op, prm, weightsCalculator, pressureIndex, comm);
520  }
521 
522  // Actually adds the creator.
523  void doAddCreator(const std::string& type, Creator c)
524  {
525  creators_[type] = c;
526  }
527 
528  // Actually adds the creator.
529  void doAddCreator(const std::string& type, ParCreator c)
530  {
531  parallel_creators_[type] = c;
532  }
533 
534  // This map contains the whole factory, i.e. all the Creators.
535  std::map<std::string, Creator> creators_;
536  std::map<std::string, ParCreator> parallel_creators_;
537 };
538 
539 } // namespace Dune
540 
541 #endif // OPM_PRECONDITIONERFACTORY_HEADER
The AMG preconditioner.
Definition: amgcpr.hh:67
Definition: PreconditionerWithUpdate.hpp:40
A two-step version of an overlapping Schwarz preconditioner using one step ILU0 as.
Definition: ParallelOverlappingILU0.hpp:611
This is an object factory for creating preconditioners.
Definition: PreconditionerFactory.hpp:51
static PrecPtr create(const Operator &op, const PropertyTree &prm, const std::function< Vector()> &weightsCalculator, const Comm &comm, std::size_t pressureIndex=std::numeric_limits< std::size_t >::max())
Create a new parallel preconditioner and return a pointer to it.
Definition: PreconditionerFactory.hpp:84
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:62
static void addCreator(const std::string &type, ParCreator creator)
Add a creator for a parallel preconditioner to the PreconditionerFactory.
Definition: PreconditionerFactory.hpp:120
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.hpp:71
std::shared_ptr< Dune::PreconditionerWithUpdate< Vector, Vector > > PrecPtr
The type of pointer returned by create().
Definition: PreconditionerFactory.hpp:58
static PrecPtr create(const Operator &op, const PropertyTree &prm, const Comm &comm, std::size_t pressureIndex=std::numeric_limits< std::size_t >::max())
Create a new parallel preconditioner and return a pointer to it.
Definition: PreconditionerFactory.hpp:96
static void addCreator(const std::string &type, Creator creator)
Add a creator for a serial preconditioner to the PreconditionerFactory.
Definition: PreconditionerFactory.hpp:108
typename Operator::matrix_type Matrix
Linear algebra types.
Definition: PreconditionerFactory.hpp:54
Definition: PropertyTree.hpp:37
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27
MILU_VARIANT
Definition: ParallelOverlappingILU0.hpp:47
@ ILU
Do not perform modified ILU.