MueLu  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
MueLu_ShiftedLaplacian_decl.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // MueLu: A package for multigrid based preconditioning
4 //
5 // Copyright 2012 NTESS and the MueLu contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
10 #ifndef MUELU_SHIFTEDLAPLACIAN_DECL_HPP
11 #define MUELU_SHIFTEDLAPLACIAN_DECL_HPP
12 
13 // Xpetra
14 #include <Xpetra_Matrix_fwd.hpp>
17 #include <Xpetra_TpetraMultiVector.hpp>
18 
19 // MueLu
20 #include "MueLu.hpp"
21 #include "MueLu_ConfigDefs.hpp"
22 
23 #if defined(HAVE_MUELU_IFPACK2)
24 
25 #include <MueLu_BaseClass.hpp>
32 #include <MueLu_Hierarchy_fwd.hpp>
34 #include <MueLu_PFactory_fwd.hpp>
35 #include <MueLu_PgPFactory_fwd.hpp>
36 #include <MueLu_RAPFactory_fwd.hpp>
38 #include <MueLu_SaPFactory_fwd.hpp>
40 #include <MueLu_ShiftedLaplacianOperator.hpp>
46 #include <MueLu_Utilities_fwd.hpp>
47 
48 // Belos
49 #include <BelosConfigDefs.hpp>
50 #include <BelosLinearProblem.hpp>
51 #include <BelosSolverFactory.hpp>
52 #include <BelosTpetraAdapter.hpp>
53 
54 #include "Kokkos_Core.hpp"
55 
56 namespace MueLu {
57 
67 template <class Scalar = DefaultScalar,
70  class Node = DefaultNode>
71 class ShiftedLaplacian : public BaseClass {
72 #undef MUELU_SHIFTEDLAPLACIAN_SHORT
73 #include "MueLu_UseShortNames.hpp"
74 
78  typedef Belos::LinearProblem<SC, TMV, OP> LinearProblem;
79  typedef Belos::SolverManager<SC, TMV, OP> SolverManager;
80  typedef Belos::SolverFactory<SC, TMV, OP> SolverFactory;
81 
82  public:
83  /*
84  FIXME 26-June-2015 JJH: This contructor is setting numerous defaults. However, they don't match the defaults
85  FIXME int the method setParameters(). There also isn't any parameter validation that I can see.
86  */
87 
90  : numPDEs_(1)
91  , Smoother_("schwarz")
92  , Aggregation_("uncoupled")
93  , Nullspace_("constant")
94  , numLevels_(5)
95  , coarseGridSize_(100)
96  , omega_(2.0 * Kokkos::numbers::pi_v<double>)
97  , iters_(500)
98  , blksize_(1)
99  , tol_(1.0e-4)
100  , nsweeps_(5)
101  , ncycles_(1)
102  , cycles_(8)
103  , subiters_(10)
104  , option_(1)
105  , nproblems_(0)
106  , solverType_(1)
107  , restart_size_(100)
108  , recycle_size_(25)
109  , smoother_sweeps_(4)
110  , smoother_damping_((SC)1.0)
111  , krylov_type_(1)
112  , krylov_iterations_(5)
114  , ilu_leveloffill_(5.0)
115  , ilu_abs_thresh_(0.0)
116  , ilu_rel_thresh_(1.0)
117  , ilu_diagpivotthresh_(0.1)
118  , ilu_drop_tol_(0.01)
119  , ilu_fill_tol_(0.01)
120  , ilu_relax_val_(1.0)
121  , ilu_rowperm_("LargeDiag")
122  , ilu_colperm_("COLAMD")
123  , ilu_drop_rule_("DROP_BASIC")
124  , ilu_normtype_("INF_NORM")
125  , ilu_milutype_("SILU")
126  , schwarz_overlap_(0)
127  , schwarz_usereorder_(true)
128  , schwarz_combinemode_(Tpetra::ADD)
129  , schwarz_ordermethod_("rcm")
130  , GridTransfersExist_(false)
131  , isSymmetric_(true) {}
132 
133  // Destructor
134  virtual ~ShiftedLaplacian();
135 
136  // Parameters
138 
139  // Set matrices
140  void setProblemMatrix(RCP<Matrix>& A);
144  void setstiff(RCP<Matrix>& K);
146  void setmass(RCP<Matrix>& M);
148  void setcoords(RCP<MultiVector>& Coords);
149  void setNullSpace(RCP<MultiVector> NullSpace);
150  void setLevelShifts(std::vector<Scalar> levelshifts);
151 
152  // initialize: set parameters and factories, construct
153  // prolongation and restriction matrices
154  void initialize();
155  // setupFastRAP: setup hierarchy with
156  // prolongators of the stiffness matrix
157  // constant complex shifts
158  void setupFastRAP();
159  // setupSlowRAP: setup hierarchy with
160  // prolongators of the stiffness matrix
161  // variable complex shifts
162  void setupSlowRAP();
163  // setupNormalRAP: setup hierarchy with
164  // prolongators of the preconditioning matrix
165  void setupNormalRAP();
166  // setupSolver: initialize Belos solver
167  void setupSolver();
168  // resetLinearProblem: for multiple frequencies;
169  // reset the Belos operator if the frequency changes
170  void resetLinearProblem();
171 
172  // Solve phase
173  Belos::ReturnType solve(const RCP<TMV> B, RCP<TMV>& X);
174  void multigrid_apply(const RCP<MultiVector> B,
175  RCP<MultiVector>& X);
178  int GetIterations();
180 
182 
183  private:
184  // Problem options
185  // numPDEs_ -> number of DOFs at each node
186  int numPDEs_;
187 
188  // Multigrid options
189  // numLevels_ -> number of Multigrid levels
190  // coarseGridSize_ -> size of coarsest grid (if current level has less DOFs, stop coarsening)
193 
194  // Shifted Laplacian/Helmholtz parameters
195  double omega_;
196  std::vector<SC> levelshifts_;
197 
198  // Krylov solver inputs
199  // iters -> max number of iterations
200  // tol -> residual tolerance
202  double tol_;
206 
207  // Smoother parameters
219  std::string schwarz_ordermethod_;
220 
221  // flags for setup
224 
225  // Xpetra matrices
226  // K_ -> stiffness matrix
227  // M_ -> mass matrix
228  // A_ -> Problem matrix
229  // P_ -> Preconditioning matrix
232 
233  // Multigrid Hierarchy
235 
236  // Factories and prototypes
251  std::string precType_;
253 
254  // Operator and Preconditioner
257 
258  // Belos Linear Problem and Solver
263 };
264 
265 } // namespace MueLu
266 
267 #define MUELU_SHIFTEDLAPLACIAN_SHORT
268 
269 #endif // if defined(HAVE_MUELU_IFPACK2) and defined(HAVE_MUELU_TPETRA)
270 
271 #endif // MUELU_SHIFTEDLAPLACIAN_DECL_HPP
RCP< SmootherFactory > coarsestSmooFact_
MueLu::DefaultLocalOrdinal LocalOrdinal
Tpetra::MultiVector< SC, LO, GO, NO > TMV
Tpetra::KokkosClassic::DefaultNode::DefaultNodeType DefaultNode
void multigrid_apply(const RCP< MultiVector > B, RCP< MultiVector > &X)
Belos::ReturnType solve(const RCP< TMV > B, RCP< TMV > &X)
Belos::SolverFactory< SC, TMV, OP > SolverFactory
RCP< CoarseMapFactory > CoarseMapfact_
RCP< CoalesceDropFactory > Dropfact_
MueLu::DefaultNode Node
void setParameters(Teuchos::RCP< Teuchos::ParameterList > paramList)
Tpetra::Operator< SC, LO, GO, NO > OP
Teuchos::ScalarTraits< Scalar >::magnitudeType GetResidual()
Belos::SolverManager< SC, TMV, OP > SolverManager
RCP< SmootherPrototype > smooProto_
MueLu::DefaultScalar Scalar
Tpetra::Details::DefaultTypes::scalar_type DefaultScalar
MueLu::DefaultGlobalOrdinal GlobalOrdinal
void setLevelShifts(std::vector< Scalar > levelshifts)
Tpetra::Vector< SC, LO, GO, NO > TVEC
RCP< Teuchos::ParameterList > BelosList_
Teuchos::ParameterList coarsestSmooList_
void setcoords(RCP< MultiVector > &Coords)
RCP< UncoupledAggregationFactory > UCaggfact_
Base class for MueLu classes.
RCP< SmootherPrototype > coarsestSmooProto_
void setProblemMatrix(RCP< Matrix > &A)
Scalar SC
RCP< Tpetra::CrsMatrix< SC, LO, GO, NO > > TpetraA_
RCP< TentativePFactory > TentPfact_
Shifted Laplacian Helmholtz solver.
RCP< MueLu::ShiftedLaplacianOperator< SC, LO, GO, NO > > MueLuOp_
Belos::LinearProblem< SC, TMV, OP > LinearProblem
RCP< AmalgamationFactory > Amalgfact_
void setNullSpace(RCP< MultiVector > NullSpace)
void setPreconditioningMatrix(RCP< Matrix > &P)