46#ifndef MUELU_MULTIVECTORTRANSFER_FACTORY_DEF_HPP
47#define MUELU_MULTIVECTORTRANSFER_FACTORY_DEF_HPP
50#include "Xpetra_MultiVectorFactory.hpp"
57 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
59 RCP<ParameterList> validParamList = rcp(
new ParameterList());
61 validParamList->set< std::string > (
"Vector name",
"undefined",
"Name of the vector that will be transferred on the coarse grid (level key)");
62 validParamList->set< RCP<const FactoryBase> >(
"Vector factory", Teuchos::null,
"Factory of the vector");
63 validParamList->set< RCP<const FactoryBase> >(
"R", Teuchos::null,
"Factory of the transfer operator (restriction)");
65 return validParamList;
68 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
70 SetParameter(
"Vector name", ParameterEntry(vectorName));
73 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
75 const ParameterList & pL = GetParameterList();
76 std::string vectorName = pL.get<std::string>(
"Vector name");
78 fineLevel.
DeclareInput(vectorName, GetFactory(
"Vector factory").get(),
this);
79 Input(coarseLevel,
"R");
82 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
86 const ParameterList & pL = GetParameterList();
87 std::string vectorName = pL.get<std::string>(
"Vector name");
89 RCP<MultiVector> fineVector = fineLevel.
Get< RCP<MultiVector> >(vectorName, GetFactory(
"Vector factory").get());
90 RCP<Matrix> transferOp = Get<RCP<Matrix> >(coarseLevel,
"R");
92 RCP<MultiVector> coarseVector = MultiVectorFactory::Build(transferOp->getRangeMap(), fineVector->getNumVectors());
93 GetOStream(
Runtime0) <<
"Transferring multivector \"" << vectorName <<
"\"" << std::endl;
95 RCP<MultiVector> onesVector = MultiVectorFactory::Build(transferOp->getDomainMap(), 1);
96 onesVector->putScalar(Teuchos::ScalarTraits<Scalar>::one());
97 RCP<MultiVector> rowSumVector = MultiVectorFactory::Build(transferOp->getRangeMap(), 1);
98 transferOp->apply(*onesVector, *rowSumVector);
99 transferOp->apply(*fineVector, *coarseVector);
101 if (vectorName ==
"Coordinates")
102 TEUCHOS_TEST_FOR_EXCEPTION(
true,
Exceptions::RuntimeError,
"Use CoordinatesTransferFactory to transfer coordinates instead of MultiVectorTransferFactory.");
104 Set<RCP<MultiVector> >(coarseLevel, vectorName, coarseVector);
108 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
113 ArrayRCP<SC> expandCoord(coordinates.size()*blksize);
115 for(
int i=0; i<coordinates.size(); i++) {
116 for(
int j=0; j< blksize; j++) {
117 expandCoord[i*blksize + j] = coordinates[i];
MueLu::DefaultLocalOrdinal LocalOrdinal
Exception throws to report errors in the internal logical of the program.
Timer to be used in factories. Similar to Monitor but with additional timers.
Class that holds all level-specific information.
void DeclareInput(const std::string &ename, const FactoryBase *factory, const FactoryBase *requestedBy=NoFactory::get())
Callback from FactoryBase::CallDeclareInput() and FactoryBase::DeclareInput()
T & Get(const std::string &ename, const FactoryBase *factory=NoFactory::get())
Get data without decrementing associated storage counter (i.e., read-only access)....
void DeclareInput(Level &finelevel, Level &coarseLevel) const
Specifies the data that this class needs, and the factories that generate that data.
void Build(Level &fineLevel, Level &coarseLevel) const
Build an object with this factory.
MultiVectorTransferFactory()
Constructor.
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.
static ArrayRCP< SC > expandCoordinates(ArrayRCP< SC > coord, LocalOrdinal blksize)
Namespace for MueLu classes and methods.
@ Runtime0
One-liner description of what is happening.