87 RCP<Matrix> A = Get< RCP<Matrix> > (fineLevel,
"A");
88 RCP<MultiVector> fineNullspace = Get< RCP<MultiVector> > (fineLevel,
"Nullspace");
90 TEUCHOS_TEST_FOR_EXCEPTION(A->GetFixedBlockSize() != 1,
Exceptions::RuntimeError,
"Block size > 1 has not been implemented");
92 const Teuchos::ParameterList& pL = GetParameterList();
94 std::string mapFile = pL.get<std::string>(
"mapFileName");
95 RCP<const Map> rowMap = A->getRowMap();
96 RCP<const Map> coarseMap = Xpetra::IO<Scalar, LocalOrdinal, GlobalOrdinal, Node>::ReadMap(mapFile, rowMap->lib(), rowMap->getComm());
97 Set(coarseLevel,
"CoarseMap", coarseMap);
99 std::string matrixFile = pL.get<std::string>(
"matrixFileName");
100 RCP<Matrix> P = Xpetra::IO<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Read(matrixFile, rowMap, coarseMap, coarseMap, rowMap);
102 Set(coarseLevel,
"P", P);
105 RCP<Matrix> P1 = Xpetra::MatrixMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Multiply(*A,
false, *P,
false);
106 P = Xpetra::IO<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Read(matrixFile, rowMap, P1->getColMap(), coarseMap, rowMap);
107 Set(coarseLevel,
"P", P);
110 RCP<MultiVector> coarseNullspace = MultiVectorFactory::Build(coarseMap, fineNullspace->getNumVectors());
111 P->apply(*fineNullspace, *coarseNullspace, Teuchos::TRANS, Teuchos::ScalarTraits<SC>::one(), Teuchos::ScalarTraits<SC>::zero());
112 Set(coarseLevel,
"Nullspace", coarseNullspace);
115 size_t n = Teuchos::as<size_t>(sqrt(coarseMap->getGlobalNumElements()));
116 TEUCHOS_TEST_FOR_EXCEPTION(n*n != coarseMap->getGlobalNumElements(),
Exceptions::RuntimeError,
"Unfortunately, this is not the case, don't know what to do");
118 RCP<MultiVector> coarseCoords = MultiVectorFactory::Build(coarseMap, 2);
119 ArrayRCP<Scalar> x = coarseCoords->getDataNonConst(0), y = coarseCoords->getDataNonConst(1);
120 for (
size_t LID = 0; LID < coarseMap->getLocalNumElements(); ++LID) {
121 GlobalOrdinal GID = coarseMap->getGlobalElement(LID) - coarseMap->getIndexBase();
126 Set(coarseLevel,
"Coordinates", coarseCoords);
129 RCP<ParameterList> params = rcp(
new ParameterList());
130 params->set(
"printLoadBalancingInfo",
true);