46#ifndef MUELU_MAXWELL1_DECL_HPP
47#define MUELU_MAXWELL1_DECL_HPP
56#include "MueLu_TrilinosSmoother.hpp"
66#include "Xpetra_Map_fwd.hpp"
67#include "Xpetra_Matrix_fwd.hpp"
68#include "Xpetra_MatrixFactory_fwd.hpp"
69#include "Xpetra_MultiVectorFactory_fwd.hpp"
70#include "Xpetra_VectorFactory_fwd.hpp"
71#include "Xpetra_CrsMatrixWrap_fwd.hpp"
87#undef MUELU_MAXWELL1_SHORT
92 typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType
magnitudeType;
93 typedef typename Teuchos::ScalarTraits<Scalar>::coordinateType
coordinateType;
106 Maxwell1(Teuchos::RCP<Hierarchy> H11, Teuchos::RCP<Hierarchy> H22) :
124 const Teuchos::RCP<Matrix> & D0_Matrix,
125 const Teuchos::RCP<MultiVector> & Nullspace,
126 const Teuchos::RCP<RealValuedMultiVector> & Coords,
127 Teuchos::ParameterList& List,
128 bool ComputePrec =
true):
131 RCP<Matrix> Kn_Matrix;
132 initialize(D0_Matrix,Kn_Matrix,Nullspace,Coords,List);
146 const Teuchos::RCP<Matrix> & D0_Matrix,
147 const Teuchos::RCP<Matrix> & Kn_Matrix,
148 const Teuchos::RCP<MultiVector> & Nullspace,
149 const Teuchos::RCP<RealValuedMultiVector> & Coords,
150 Teuchos::ParameterList& List,
151 bool ComputePrec =
true):
154 initialize(D0_Matrix,Kn_Matrix,Nullspace,Coords,List);
169 const Teuchos::RCP<Matrix> & D0_Matrix,
170 const Teuchos::RCP<Matrix> & Kn_Matrix,
171 const Teuchos::RCP<MultiVector> & Nullspace,
172 const Teuchos::RCP<RealValuedMultiVector> & Coords,
173 Teuchos::ParameterList& List,
const Teuchos::RCP<Matrix> & GmhdA_Matrix,
174 bool ComputePrec =
true):
177 initialize(D0_Matrix,Kn_Matrix,Nullspace,Coords,List);
191 Teuchos::ParameterList& List,
192 bool ComputePrec =
true):
196 RCP<MultiVector> Nullspace = List.get<RCP<MultiVector> >(
"Nullspace", Teuchos::null);
197 RCP<RealValuedMultiVector> Coords = List.get<RCP<RealValuedMultiVector> >(
"Coordinates", Teuchos::null);
198 RCP<Matrix> D0_Matrix = List.get<RCP<Matrix> >(
"D0");
199 RCP<Matrix> Kn_Matrix;
200 if (List.isType<RCP<Matrix> >(
"Kn"))
201 Kn_Matrix = List.get<RCP<Matrix> >(
"Kn");
203 initialize(D0_Matrix,Kn_Matrix,Nullspace,Coords,List);
205 if (SM_Matrix != Teuchos::null)
227 void compute(
bool reuse=
false);
230 void resetMatrix(Teuchos::RCP<Matrix> SM_Matrix_new,
bool ComputePrec=
true);
235 void apply (
const MultiVector& X, MultiVector& Y,
236 Teuchos::ETransp mode = Teuchos::NO_TRANS,
237 Scalar alpha = Teuchos::ScalarTraits<Scalar>::one(),
238 Scalar beta = Teuchos::ScalarTraits<Scalar>::zero())
const;
243 void describe(Teuchos::FancyOStream &out,
const Teuchos::EVerbosityLevel verbLevel = Teuchos::VERB_HIGH)
const;
247 const MultiVector & B,
248 MultiVector & R)
const {
249 using STS = Teuchos::ScalarTraits<Scalar>;
250 R.update(STS::one(),B,STS::zero());
251 this->
apply (X, R, Teuchos::NO_TRANS, -STS::one(), STS::one());
271 void initialize(
const Teuchos::RCP<Matrix> & D0_Matrix,
272 const Teuchos::RCP<Matrix> & Kn_Matrix,
273 const Teuchos::RCP<MultiVector> & Nullspace,
274 const Teuchos::RCP<RealValuedMultiVector> & Coords,
275 Teuchos::ParameterList& List);
288 void dump(
const Matrix& A, std::string name)
const;
291 void dump(
const MultiVector& X, std::string name)
const;
297 void dump(
const Teuchos::ArrayRCP<bool>& v, std::string name)
const;
300 void dump(
const Kokkos::View<bool*, typename Node::device_type>& v, std::string name)
const;
341#define MUELU_MAXWELL1_SHORT
MueLu::DefaultLocalOrdinal LocalOrdinal
MueLu::DefaultScalar Scalar
MueLu::DefaultGlobalOrdinal GlobalOrdinal
Provides methods to build a multigrid hierarchy and apply multigrid cycles.
Preconditioner (wrapped as a Xpetra::Operator) for Maxwell's equations in curl-curl form.
Teuchos::ScalarTraits< Scalar >::coordinateType coordinateType
Xpetra::MultiVector< coordinateType, LO, GO, NO > RealValuedMultiVector
bool useKokkos_
Some options.
Teuchos::RCP< Matrix > generate_kn() const
Generates the Kn matrix.
Teuchos::RCP< Hierarchy > Hierarchy11_
Two hierarchies: one for the (1,1)-block, another for the (2,2)-block.
Teuchos::RCP< Hierarchy > HierarchyGmhd_
Teuchos::RCP< MultiVector > residualFine_
Teuchos::ArrayRCP< bool > BCcols_
virtual ~Maxwell1()
Destructor.
Teuchos::ArrayRCP< bool > BCrows_
Teuchos::RCP< Matrix > SM_Matrix_
Various matrices.
Teuchos::RCP< Hierarchy > Hierarchy22_
Maxwell1(const Teuchos::RCP< Matrix > &SM_Matrix, const Teuchos::RCP< Matrix > &D0_Matrix, const Teuchos::RCP< MultiVector > &Nullspace, const Teuchos::RCP< RealValuedMultiVector > &Coords, Teuchos::ParameterList &List, bool ComputePrec=true)
Teuchos::RCP< const Map > getDomainMap() const
Returns the Xpetra::Map object associated with the domain of this operator.
Maxwell1(const Teuchos::RCP< Matrix > &SM_Matrix, const Teuchos::RCP< Matrix > &D0_Matrix, const Teuchos::RCP< Matrix > &Kn_Matrix, const Teuchos::RCP< MultiVector > &Nullspace, const Teuchos::RCP< RealValuedMultiVector > &Coords, Teuchos::ParameterList &List, bool ComputePrec=true)
Teuchos::RCP< const Map > getRangeMap() const
Returns the Xpetra::Map object associated with the range of this operator.
void compute(bool reuse=false)
Setup the preconditioner.
void GMHDSetupHierarchy(Teuchos::ParameterList &List) const
Sets up hiearchy for GMHD matrices that include generalized Ohms law equations.
mode_type
Execution modes.
Teuchos::ParameterList parameterList_
ParameterLists.
Teuchos::RCP< MultiVector > residual22_
Teuchos::RCP< RealValuedMultiVector > Coords_
Coordinates.
Teuchos::RCP< Matrix > D0_Matrix_
Teuchos::ScalarTraits< Scalar >::magnitudeType magnitudeType
void applyInverseRefMaxwellAdditive(const MultiVector &RHS, MultiVector &X) const
apply RefMaxwell additive 2x2 style cycle
Teuchos::RCP< MultiVector > Nullspace_
Nullspace.
Kokkos::View< bool *, typename Node::device_type > BCdomainKokkos_
Teuchos::RCP< MultiVector > update11c_
Kokkos::View< bool *, typename Node::device_type > BCcolsKokkos_
Maxwell1(Teuchos::RCP< Hierarchy > H11, Teuchos::RCP< Hierarchy > H22)
Constructor with Hierarchies.
void applyInverseStandard(const MultiVector &RHS, MultiVector &X) const
apply standard Maxwell1 cycle
void initialize(const Teuchos::RCP< Matrix > &D0_Matrix, const Teuchos::RCP< Matrix > &Kn_Matrix, const Teuchos::RCP< MultiVector > &Nullspace, const Teuchos::RCP< RealValuedMultiVector > &Coords, Teuchos::ParameterList &List)
RCP< Matrix > P11_
Temporary memory (cached vectors for RefMaxwell-style)
void dumpCoords(const RealValuedMultiVector &X, std::string name) const
dump out real-valued multivector
Teuchos::RCP< MultiVector > update22_
Teuchos::ArrayRCP< bool > BCdomain_
Teuchos::ParameterList precList11_
Maxwell1(const Teuchos::RCP< Matrix > &SM_Matrix, const Teuchos::RCP< Matrix > &D0_Matrix, const Teuchos::RCP< Matrix > &Kn_Matrix, const Teuchos::RCP< MultiVector > &Nullspace, const Teuchos::RCP< RealValuedMultiVector > &Coords, Teuchos::ParameterList &List, const Teuchos::RCP< Matrix > &GmhdA_Matrix, bool ComputePrec=true)
void residual(const MultiVector &X, const MultiVector &B, MultiVector &R) const
Compute a residual R = B - (*this) * X.
const Teuchos::RCP< Matrix > & getJacobian() const
Returns Jacobian matrix SM.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::VERB_HIGH) const
Teuchos::RCP< Matrix > GmhdA_Matrix_
void dump(const Matrix &A, std::string name) const
dump out matrix
Maxwell1(const Teuchos::RCP< Matrix > &SM_Matrix, Teuchos::ParameterList &List, bool ComputePrec=true)
Teuchos::ParameterList precList22_
void resetMatrix(Teuchos::RCP< Matrix > SM_Matrix_new, bool ComputePrec=true)
Reset system matrix.
Teuchos::RCP< Matrix > Kn_Matrix_
Teuchos::RCP< MultiVector > residual11c_
bool hasTransposeApply() const
Indicates whether this operator supports applying the adjoint operator.
Kokkos::View< bool *, typename Node::device_type > BCrowsKokkos_
Vectors for BCs.
Teuchos::RCP< Teuchos::TimeMonitor > getTimer(std::string name, RCP< const Teuchos::Comm< int > > comm=Teuchos::null) const
get a (synced) timer
void setParameters(Teuchos::ParameterList &list)
Set parameters.
void apply(const MultiVector &X, MultiVector &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=Teuchos::ScalarTraits< Scalar >::one(), Scalar beta=Teuchos::ScalarTraits< Scalar >::zero()) const
void allocateMemory(int numVectors) const
allocate multivectors for solve
Verbose class for MueLu classes.
Namespace for MueLu classes and methods.