AMDiS  0.3
The Adaptive Multi-Dimensional Simulation Toolbox
Creators.hpp
1 #pragma once
2 
3 #include <memory>
4 #include <utility>
5 
6 #include <dune/common/typeutilities.hh>
7 #include <dune/istl/operators.hh>
8 #include <dune/istl/scalarproducts.hh>
9 #include <dune/istl/solvercategory.hh>
10 
11 #if HAVE_MPI
12 #include <dune/istl/novlpschwarz.hh>
13 #include <dune/istl/schwarz.hh>
14 #endif
15 
16 #include <amdis/Output.hpp>
17 #include <amdis/linearalgebra/istl/PreconWrapper.hpp>
18 
19 namespace AMDiS
20 {
22  template <class X>
24  {
25  public:
26  using Interface = Dune::ScalarProduct<X>;
27 
28  template <class Comm>
29  static std::unique_ptr<Interface> create(Dune::SolverCategory::Category cat, Comm const& comm)
30  {
31  return create_impl(cat, comm.get(), Dune::PriorityTag<10>{});
32  }
33 
34  private:
35  static std::unique_ptr<Interface>
36  create_impl(Dune::SolverCategory::Category cat, Dune::Amg::SequentialInformation, Dune::PriorityTag<2>)
37  {
38  assert(cat == Dune::SolverCategory::sequential);
39  return std::make_unique< Dune::SeqScalarProduct<X> >();
40  }
41 
42  template <class Comm>
43  static std::unique_ptr<Interface>
44  create_impl(Dune::SolverCategory::Category cat, Comm const& comm, Dune::PriorityTag<1>)
45  {
46  switch (cat)
47  {
48  case Dune::SolverCategory::sequential:
49  return std::make_unique< Dune::SeqScalarProduct<X> >();
50 #if HAVE_MPI
51  case Dune::SolverCategory::overlapping:
52  return std::make_unique< Dune::OverlappingSchwarzScalarProduct<X,Comm> >(comm);
53  case Dune::SolverCategory::nonoverlapping:
54  return std::make_unique< Dune::NonoverlappingSchwarzScalarProduct<X,Comm> >(comm);
55 #endif
56  default:
57  error_exit("Invalid solver category for scalar product\n");
58  return nullptr; // avoid warnings
59  }
60  }
61  };
62 
64  template <class M, class X, class Y>
66  {
67  public:
68  using Interface = Dune::AssembledLinearOperator<M,X,Y>;
69 
70  template <class Comm>
71  static std::unique_ptr<Interface> create(Dune::SolverCategory::Category cat, M const& mat, Comm const& comm)
72  {
73  return create_impl(cat, mat, comm.get(), Dune::PriorityTag<10>{});
74  }
75 
76  private:
77  static std::unique_ptr<Interface>
78  create_impl(Dune::SolverCategory::Category cat, M const& mat, Dune::Amg::SequentialInformation, Dune::PriorityTag<2>)
79  {
80  assert(cat == Dune::SolverCategory::sequential);
81  return std::make_unique< Dune::MatrixAdapter<M,X,Y> >(mat);
82  }
83 
84  template <class Comm>
85  static std::unique_ptr<Interface>
86  create_impl(Dune::SolverCategory::Category cat, M const& mat, Comm const& comm, Dune::PriorityTag<1>)
87  {
88  switch (cat)
89  {
90  case Dune::SolverCategory::sequential:
91  return std::make_unique< Dune::MatrixAdapter<M,X,Y> >(mat);
92 #if HAVE_MPI
93  case Dune::SolverCategory::overlapping:
94  return std::make_unique< Dune::OverlappingSchwarzOperator<M,X,Y,Comm> >(mat, comm);
95  case Dune::SolverCategory::nonoverlapping:
96  return std::make_unique< Dune::NonoverlappingSchwarzOperator<M,X,Y,Comm> >(mat, comm);
97 #endif
98  default:
99  error_exit("Invalid solver category for linear operator\n");
100  return nullptr; // avoid warnings
101  }
102  }
103  };
104 
105 
107  template <class X, class Y>
109  {
110  public:
111  using Interface = Dune::Preconditioner<X,Y>;
112 
113  template <class Comm>
114  static std::unique_ptr<Interface> create(Dune::SolverCategory::Category cat, std::unique_ptr<Interface> prec, Comm const& comm)
115  {
116  return create_impl(cat, std::move(prec), comm.get(), Dune::PriorityTag<10>{});
117  }
118 
119  private:
120  static std::unique_ptr<Interface>
121  create_impl(Dune::SolverCategory::Category cat, std::unique_ptr<Interface> prec, Dune::Amg::SequentialInformation, Dune::PriorityTag<2>)
122  {
123  assert(cat == Dune::SolverCategory::sequential);
124  return std::move(prec);
125  }
126 
127  template <class Comm>
128  static std::unique_ptr<Interface>
129  create_impl(Dune::SolverCategory::Category cat, std::unique_ptr<Interface> prec, Comm const& comm, Dune::PriorityTag<1>)
130  {
131  switch (cat)
132  {
133  case Dune::SolverCategory::sequential:
134  return std::move(prec);
135 #if HAVE_MPI
136  case Dune::SolverCategory::overlapping:
137  {
138  using BP = Dune::BlockPreconditioner<X,Y,Comm>;
139  return std::make_unique< PreconWrapper<BP,Interface> >(std::move(prec), comm);
140  }
141  case Dune::SolverCategory::nonoverlapping:
142  {
143  using BP = Dune::NonoverlappingBlockPreconditioner<Comm,Interface>;
144  return std::make_unique< PreconWrapper<BP,Interface> >(std::move(prec), comm);
145  }
146 #endif
147  default:
148  error_exit("Invalid solver category for parallel preconditioner\n");
149  return nullptr; // avoid warnings
150  }
151  }
152  };
153 
154 } // end namespace AMDiS
void error_exit(std::string const &str, Args &&... args)
print a message and exit
Definition: Output.hpp:142
Contains all classes needed for solving linear and non linear equation systems.
Definition: AdaptBase.hpp:6
auto cat(Dune::TypeTree::HybridTreePath< S... > const &tp0, Dune::TypeTree::HybridTreePath< T... > const &tp1)
Concatenate two treepaths.
Definition: TreePath.hpp:209
Creator to create ScalarProduct objects.
Definition: Creators.hpp:23
Creator to create Linear Operator objects.
Definition: Creators.hpp:65
Creator to create parallel Preconditioners.
Definition: Creators.hpp:108