46 #ifndef MUELU_SHIFTEDLAPLACIAN_DEF_HPP 47 #define MUELU_SHIFTEDLAPLACIAN_DEF_HPP 51 #if defined(HAVE_MUELU_IFPACK2) and defined(HAVE_MUELU_TPETRA) 53 #include <MueLu_AmalgamationFactory.hpp> 54 #include <MueLu_CoalesceDropFactory.hpp> 55 #include <MueLu_CoarseMapFactory.hpp> 56 #include <MueLu_CoupledAggregationFactory.hpp> 57 #include <MueLu_CoupledRBMFactory.hpp> 58 #include <MueLu_DirectSolver.hpp> 59 #include <MueLu_GenericRFactory.hpp> 60 #include <MueLu_Hierarchy.hpp> 61 #include <MueLu_Ifpack2Smoother.hpp> 62 #include <MueLu_PFactory.hpp> 63 #include <MueLu_PgPFactory.hpp> 64 #include <MueLu_RAPFactory.hpp> 65 #include <MueLu_RAPShiftFactory.hpp> 66 #include <MueLu_SaPFactory.hpp> 67 #include <MueLu_ShiftedLaplacian.hpp> 68 #include <MueLu_ShiftedLaplacianOperator.hpp> 69 #include <MueLu_SmootherFactory.hpp> 70 #include <MueLu_SmootherPrototype.hpp> 71 #include <MueLu_TentativePFactory.hpp> 72 #include <MueLu_TransPFactory.hpp> 73 #include <MueLu_UncoupledAggregationFactory.hpp> 74 #include <MueLu_Utilities.hpp> 79 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
83 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
87 coarseGridSize_ = paramList->get(
"MueLu: coarse size", 1000);
88 numLevels_ = paramList->get(
"MueLu: levels", 3);
89 int stype = paramList->get(
"MueLu: smoother", 8);
90 if(stype==1) { Smoother_=
"jacobi"; }
91 else if(stype==2) { Smoother_=
"gauss-seidel"; }
92 else if(stype==3) { Smoother_=
"symmetric gauss-seidel"; }
93 else if(stype==4) { Smoother_=
"chebyshev"; }
94 else if(stype==5) { Smoother_=
"krylov"; }
95 else if(stype==6) { Smoother_=
"ilut"; }
96 else if(stype==7) { Smoother_=
"riluk"; }
97 else if(stype==8) { Smoother_=
"schwarz"; }
98 else if(stype==9) { Smoother_=
"superilu"; }
99 else if(stype==10) { Smoother_=
"superlu"; }
100 else { Smoother_=
"schwarz"; }
101 smoother_sweeps_ = paramList->get(
"MueLu: sweeps", 5);
102 smoother_damping_ = paramList->get(
"MueLu: relax val", 1.0);
103 ncycles_ = paramList->get(
"MueLu: cycles", 1);
104 iters_ = paramList->get(
"MueLu: iterations", 500);
105 solverType_ = paramList->get(
"MueLu: solver type", 1);
106 restart_size_ = paramList->get(
"MueLu: restart size", 100);
107 recycle_size_ = paramList->get(
"MueLu: recycle size", 25);
108 isSymmetric_ = paramList->get(
"MueLu: symmetric",
true);
109 ilu_leveloffill_ = paramList->get(
"MueLu: level-of-fill", 5);
110 ilu_abs_thresh_ = paramList->get(
"MueLu: abs thresh", 0.0);
111 ilu_rel_thresh_ = paramList->get(
"MueLu: rel thresh", 1.0);
112 ilu_diagpivotthresh_ = paramList->get(
"MueLu: piv thresh", 0.1);
113 ilu_drop_tol_ = paramList->get(
"MueLu: drop tol", 0.01);
114 ilu_fill_tol_ = paramList->get(
"MueLu: fill tol", 0.01);
115 schwarz_overlap_ = paramList->get(
"MueLu: overlap", 0);
116 schwarz_usereorder_ = paramList->get(
"MueLu: use reorder",
true);
117 int combinemode = paramList->get(
"MueLu: combine mode", 1);
118 if(combinemode==0) { schwarz_combinemode_ = Tpetra::ZERO; }
119 else { schwarz_combinemode_ = Tpetra::ADD; }
120 tol_ = paramList->get(
"MueLu: tolerance", 0.001);
124 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
128 if(A_!=Teuchos::null)
130 #ifdef HAVE_MUELU_TPETRA_INST_INT_INT 131 if(LinearProblem_!=Teuchos::null)
132 LinearProblem_ -> setOperator ( TpetraA_ );
134 TEUCHOS_TEST_FOR_EXCEPTION(
true,
Exceptions::RuntimeError,
"ShiftedLaplacian only available with Tpetra and GO=int enabled.");
139 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
143 #ifdef HAVE_MUELU_TPETRA_INST_INT_INT 144 if(LinearProblem_!=Teuchos::null)
145 LinearProblem_ -> setOperator ( TpetraA_ );
150 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
154 GridTransfersExist_=
false;
158 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
161 RCP< Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > Atmp
162 = rcp(
new Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(TpetraP) );
163 P_= rcp(
new Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>(Atmp) );
164 GridTransfersExist_=
false;
168 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
175 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
178 RCP< Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > Atmp
179 = rcp(
new Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(TpetraK) );
180 K_= rcp(
new Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>(Atmp) );
184 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
191 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
194 RCP< Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > Atmp
195 = rcp(
new Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(TpetraM) );
196 M_= rcp(
new Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>(Atmp) );
200 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
207 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
210 NullSpace_=NullSpace;
214 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
217 levelshifts_=levelshifts;
218 numLevels_=levelshifts_.size();
223 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
239 Manager_ -> SetFactory(
"UnAmalgamationInfo", Amalgfact_);
240 Teuchos::ParameterList params;
241 params.set(
"lightweight wrap",
true);
242 params.set(
"aggregation: drop scheme",
"classical");
243 Dropfact_ -> SetParameterList(params);
244 Manager_ -> SetFactory(
"Graph", Dropfact_);
245 if(Aggregation_==
"coupled") {
246 Manager_ -> SetFactory(
"Aggregates", Aggfact_ );
249 Manager_ -> SetFactory(
"Aggregates", UCaggfact_ );
251 Manager_ -> SetFactory(
"CoarseMap", CoarseMapfact_);
252 Manager_ -> SetFactory(
"Ptent", TentPfact_);
253 if(isSymmetric_==
true) {
254 Manager_ -> SetFactory(
"P", Pfact_);
255 Manager_ -> SetFactory(
"R", TransPfact_);
258 Manager_ -> SetFactory(
"P", PgPfact_);
259 Manager_ -> SetFactory(
"R", Rfact_);
264 if(Smoother_==
"jacobi") {
265 precType_ =
"RELAXATION";
266 precList_.set(
"relaxation: type",
"Jacobi");
267 precList_.set(
"relaxation: sweeps", smoother_sweeps_);
268 precList_.set(
"relaxation: damping factor", smoother_damping_);
270 else if(Smoother_==
"gauss-seidel") {
271 precType_ =
"RELAXATION";
272 precList_.set(
"relaxation: type",
"Gauss-Seidel");
273 precList_.set(
"relaxation: sweeps", smoother_sweeps_);
274 precList_.set(
"relaxation: damping factor", smoother_damping_);
276 else if(Smoother_==
"symmetric gauss-seidel") {
277 precType_ =
"RELAXATION";
278 precList_.set(
"relaxation: type",
"Symmetric Gauss-Seidel");
279 precList_.set(
"relaxation: sweeps", smoother_sweeps_);
280 precList_.set(
"relaxation: damping factor", smoother_damping_);
282 else if(Smoother_==
"chebyshev") {
283 precType_ =
"CHEBYSHEV";
285 else if(Smoother_==
"krylov") {
286 precType_ =
"KRYLOV";
287 precList_.set(
"krylov: iteration type", krylov_type_);
288 precList_.set(
"krylov: number of iterations", krylov_iterations_);
289 precList_.set(
"krylov: residual tolerance",1.0e-8);
290 precList_.set(
"krylov: block size",1);
291 precList_.set(
"krylov: preconditioner type", krylov_preconditioner_);
292 precList_.set(
"relaxation: sweeps",1);
295 else if(Smoother_==
"ilut") {
297 precList_.set(
"fact: ilut level-of-fill", ilu_leveloffill_);
298 precList_.set(
"fact: absolute threshold", ilu_abs_thresh_);
299 precList_.set(
"fact: relative threshold", ilu_rel_thresh_);
300 precList_.set(
"fact: drop tolerance", ilu_drop_tol_);
301 precList_.set(
"fact: relax value", ilu_relax_val_);
303 else if(Smoother_==
"riluk") {
305 precList_.set(
"fact: iluk level-of-fill", ilu_leveloffill_);
306 precList_.set(
"fact: absolute threshold", ilu_abs_thresh_);
307 precList_.set(
"fact: relative threshold", ilu_rel_thresh_);
308 precList_.set(
"fact: drop tolerance", ilu_drop_tol_);
309 precList_.set(
"fact: relax value", ilu_relax_val_);
311 else if(Smoother_==
"schwarz") {
312 precType_ =
"SCHWARZ";
313 precList_.set(
"schwarz: overlap level", schwarz_overlap_);
314 precList_.set(
"schwarz: combine mode", schwarz_combinemode_);
315 precList_.set(
"schwarz: use reordering", schwarz_usereorder_);
317 precList_.set(
"order_method",schwarz_ordermethod_);
318 precList_.sublist(
"schwarz: reordering list").set(
"order_method",schwarz_ordermethod_);
319 precList_.sublist(
"schwarz: subdomain solver parameters").set(
"fact: ilut level-of-fill", ilu_leveloffill_);
320 precList_.sublist(
"schwarz: subdomain solver parameters").set(
"fact: absolute threshold", ilu_abs_thresh_);
321 precList_.sublist(
"schwarz: subdomain solver parameters").set(
"fact: relative threshold", ilu_rel_thresh_);
322 precList_.sublist(
"schwarz: subdomain solver parameters").set(
"fact: drop tolerance", ilu_drop_tol_);
323 precList_.sublist(
"schwarz: subdomain solver parameters").set(
"fact: relax value", ilu_relax_val_);
325 else if(Smoother_==
"superilu") {
326 precType_ =
"superlu";
327 precList_.set(
"RowPerm", ilu_rowperm_);
328 precList_.set(
"ColPerm", ilu_colperm_);
329 precList_.set(
"DiagPivotThresh", ilu_diagpivotthresh_);
330 precList_.set(
"ILU_DropRule",ilu_drop_rule_);
331 precList_.set(
"ILU_DropTol",ilu_drop_tol_);
332 precList_.set(
"ILU_FillFactor",ilu_leveloffill_);
333 precList_.set(
"ILU_Norm",ilu_normtype_);
334 precList_.set(
"ILU_MILU",ilu_milutype_);
335 precList_.set(
"ILU_FillTol",ilu_fill_tol_);
336 precList_.set(
"ILU_Flag",
true);
338 else if(Smoother_==
"superlu") {
339 precType_ =
"superlu";
340 precList_.set(
"ColPerm", ilu_colperm_);
341 precList_.set(
"DiagPivotThresh", ilu_diagpivotthresh_);
343 #ifdef HAVE_MUELU_TPETRA_INST_INT_INT 347 #if defined(HAVE_MUELU_AMESOS2) and defined(HAVE_AMESOS2_SUPERLU) 348 coarsestSmooProto_ = rcp(
new DirectSolver(
"Superlu",coarsestSmooList_) );
349 #elif defined(HAVE_MUELU_AMESOS2) and defined(HAVE_AMESOS2_KLU2) 350 coarsestSmooProto_ = rcp(
new DirectSolver(
"Klu",coarsestSmooList_) );
351 #elif defined(HAVE_MUELU_AMESOS2) and defined(HAVE_AMESOS2_SUPERLUDIST) 352 coarsestSmooProto_ = rcp(
new DirectSolver(
"Superludist",coarsestSmooList_) );
356 coarsestSmooFact_ = rcp(
new SmootherFactory(coarsestSmooProto_, Teuchos::null) );
364 if(K_!=Teuchos::null) {
365 Manager_ -> SetFactory(
"Smoother", Teuchos::null);
366 Manager_ -> SetFactory(
"CoarseSolver", Teuchos::null);
368 if(NullSpace_!=Teuchos::null)
369 Hierarchy_ -> GetLevel(0) -> Set(
"Nullspace", NullSpace_);
370 if(isSymmetric_==
true) {
371 Hierarchy_ ->
Keep(
"P", Pfact_.get());
372 Hierarchy_ ->
Keep(
"R", TransPfact_.get());
373 Hierarchy_ -> SetImplicitTranspose(
true);
376 Hierarchy_ ->
Keep(
"P", PgPfact_.get());
377 Hierarchy_ ->
Keep(
"R", Rfact_.get());
379 Hierarchy_ ->
Keep(
"Ptent", TentPfact_.get());
380 Hierarchy_ -> SetMaxCoarseSize( coarseGridSize_ );
381 Hierarchy_ -> Setup(*Manager_, 0, numLevels_);
382 GridTransfersExist_=
true;
386 Manager_ -> SetFactory(
"Smoother", smooFact_);
387 Manager_ -> SetFactory(
"CoarseSolver", coarsestSmooFact_);
389 if(NullSpace_!=Teuchos::null)
390 Hierarchy_ -> GetLevel(0) -> Set(
"Nullspace", NullSpace_);
391 if(isSymmetric_==
true)
392 Hierarchy_ -> SetImplicitTranspose(
true);
393 Hierarchy_ -> SetMaxCoarseSize( coarseGridSize_ );
394 Hierarchy_ -> Setup(*Manager_, 0, numLevels_);
395 GridTransfersExist_=
true;
399 BelosList_ = rcp(
new Teuchos::ParameterList(
"GMRES") );
400 BelosList_ ->
set(
"Maximum Iterations",iters_ );
401 BelosList_ ->
set(
"Convergence Tolerance",tol_ );
403 BelosList_ ->
set(
"Output Frequency",1);
404 BelosList_ ->
set(
"Output Style",Belos::Brief);
405 BelosList_ ->
set(
"Num Blocks",restart_size_);
406 BelosList_ ->
set(
"Num Recycled Blocks",recycle_size_);
408 TEUCHOS_TEST_FOR_EXCEPTION(
true,
Exceptions::RuntimeError,
"ShiftedLaplacian only available with Tpetra and GO=int enabled.");
413 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
416 int numLevels = Hierarchy_ -> GetNumLevels();
417 Manager_ -> SetFactory(
"Smoother", smooFact_);
418 Manager_ -> SetFactory(
"CoarseSolver", coarsestSmooFact_);
419 Hierarchy_ -> GetLevel(0) -> Set(
"A", P_);
420 Hierarchy_ -> Setup(*Manager_, 0, numLevels);
426 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
429 int numLevels = Hierarchy_ -> GetNumLevels();
430 Acshift_->SetShifts(levelshifts_);
431 Manager_ -> SetFactory(
"Smoother", smooFact_);
432 Manager_ -> SetFactory(
"CoarseSolver", coarsestSmooFact_);
433 Manager_ -> SetFactory(
"A", Acshift_);
434 Manager_ -> SetFactory(
"K", Acshift_);
435 Manager_ -> SetFactory(
"M", Acshift_);
436 Hierarchy_ -> GetLevel(0) -> Set(
"A", P_);
437 Hierarchy_ -> GetLevel(0) -> Set(
"K", K_);
438 Hierarchy_ -> GetLevel(0) -> Set(
"M", M_);
439 Hierarchy_ -> Setup(*Manager_, 0, numLevels);
445 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
449 if( GridTransfersExist_ ==
false ) {
451 if(NullSpace_!=Teuchos::null)
452 Hierarchy_ -> GetLevel(0) -> Set(
"Nullspace", NullSpace_);
453 if(isSymmetric_==
true)
454 Hierarchy_ -> SetImplicitTranspose(
true);
455 Hierarchy_ -> SetMaxCoarseSize( coarseGridSize_ );
456 Hierarchy_ -> Setup(*Manager_, 0, numLevels_);
457 GridTransfersExist_=
true;
463 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
466 #ifdef HAVE_MUELU_TPETRA_INST_INT_INT 469 (Hierarchy_, A_, ncycles_, subiters_, option_, tol_) );
471 if(LinearProblem_==Teuchos::null)
472 LinearProblem_ = rcp(
new LinearProblem );
473 LinearProblem_ -> setOperator ( TpetraA_ );
474 LinearProblem_ -> setRightPrec( MueLuOp_ );
475 if(SolverManager_==Teuchos::null) {
476 std::string solverName;
477 SolverFactory_= rcp(
new SolverFactory() );
478 if(solverType_==1) { solverName=
"Block GMRES"; }
479 else if(solverType_==2) { solverName=
"Recycling GMRES"; }
480 else { solverName=
"Flexible GMRES"; }
481 SolverManager_ = SolverFactory_->create( solverName, BelosList_ );
482 SolverManager_ -> setProblem( LinearProblem_ );
485 TEUCHOS_TEST_FOR_EXCEPTION(
true,
Exceptions::RuntimeError,
"ShiftedLaplacian only available with Tpetra and GO=int enabled.");
489 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
492 #ifdef HAVE_MUELU_TPETRA_INST_INT_INT 493 LinearProblem_ -> setOperator ( TpetraA_ );
495 TEUCHOS_TEST_FOR_EXCEPTION(
true,
Exceptions::RuntimeError,
"ShiftedLaplacian only available with Tpetra and GO=int enabled.");
500 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
503 #ifdef HAVE_MUELU_TPETRA_INST_INT_INT 505 LinearProblem_ -> setProblem(X, B);
507 SolverManager_ -> solve();
509 TEUCHOS_TEST_FOR_EXCEPTION(
true,
Exceptions::RuntimeError,
"ShiftedLaplacian only available with Tpetra and GO=int enabled.");
515 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
520 Hierarchy_ -> Iterate(*B, *X, 1,
true, 0);
524 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
526 RCP<Tpetra::MultiVector<SC,LO,GO,NO> >& X)
528 Teuchos::RCP< Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > XpetraX
529 = Teuchos::rcp(
new Xpetra::TpetraMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>(X) );
530 Teuchos::RCP< Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > XpetraB
531 = Teuchos::rcp(
new Xpetra::TpetraMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>(B) );
533 Hierarchy_ -> Iterate(*XpetraB, *XpetraX, 1,
true, 0);
537 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
540 #ifdef HAVE_MUELU_TPETRA_INST_INT_INT 541 int numiters = SolverManager_ -> getNumIters();
544 TEUCHOS_TEST_FOR_EXCEPTION(
true,
Exceptions::RuntimeError,
"ShiftedLaplacian only available with Tpetra and GO=int enabled.");
550 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
551 typename Teuchos::ScalarTraits<Scalar>::magnitudeType
554 typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType MT;
555 #ifdef HAVE_MUELU_TPETRA_INST_INT_INT 556 MT residual = SolverManager_ -> achievedTol();
559 TEUCHOS_TEST_FOR_EXCEPTION(
true,
Exceptions::RuntimeError,
"ShiftedLaplacian only available with Tpetra and GO=int enabled.");
566 #define MUELU_SHIFTEDLAPLACIAN_SHORT 568 #endif //if defined(HAVE_MUELU_IFPACK2) and defined(HAVE_MUELU_TPETRA) 569 #endif // MUELU_SHIFTEDLAPLACIAN_DEF_HPP Generic Smoother Factory for generating the smoothers of the MG hierarchy.
This class specifies the default factory that should generate some data on a Level if the data does n...
Factory for generating coarse level map. Used by TentativePFactory.
Factory for building coarse grid matrices, when the matrix is of the form K+a*M. Useful when you want...
void multigrid_apply(const RCP< MultiVector > B, RCP< MultiVector > &X)
void setmass(RCP< Matrix > &M)
Namespace for MueLu classes and methods.
void resetLinearProblem()
Factory for building tentative prolongator.
void setParameters(Teuchos::RCP< Teuchos::ParameterList > paramList)
Class that encapsulates direct solvers. Autoselection of AmesosSmoother or Amesos2Smoother according ...
Factory for coarsening a graph with uncoupled aggregation.
Factory for building restriction operators using a prolongator factory.
Teuchos::ScalarTraits< Scalar >::magnitudeType GetResidual()
void setLevelShifts(std::vector< Scalar > levelshifts)
static RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraCrs(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
virtual ~ShiftedLaplacian()
Always keep data, even accross run. This flag is set by Level::Keep(). This flag is propagated to coa...
AmalgamationFactory for subblocks of strided map based amalgamation data.
void setcoords(RCP< MultiVector > &Coords)
Wraps an existing MueLu::Hierarchy as a Tpetra::Operator, with an optional two-level correction...
int solve(const RCP< TMV > B, RCP< TMV > &X)
Factory for creating a graph based on a given matrix.
void setProblemMatrix(RCP< Matrix > &A)
void setstiff(RCP< Matrix > &K)
Factory for building restriction operators.
Exception throws to report errors in the internal logical of the program.
Print all warning messages.
Factory for building coarse matrices.
Factory for building Petrov-Galerkin Smoothed Aggregation prolongators.
Class that encapsulates Ifpack2 smoothers.
Factory for building Smoothed Aggregation prolongators.
Factory for building uncoupled aggregates.
void setNullSpace(RCP< MultiVector > NullSpace)
Provides methods to build a multigrid hierarchy and apply multigrid cycles.
void setPreconditioningMatrix(RCP< Matrix > &P)