46 #ifndef MUELU_TENTATIVEPFACTORY_DEF_HPP 47 #define MUELU_TENTATIVEPFACTORY_DEF_HPP 64 #include "MueLu_Aggregates.hpp" 65 #include "MueLu_AmalgamationFactory.hpp" 66 #include "MueLu_AmalgamationInfo.hpp" 67 #include "MueLu_CoarseMapFactory.hpp" 68 #include "MueLu_NullspaceFactory.hpp" 69 #include "MueLu_PerfUtils.hpp" 71 #include "MueLu_Utilities.hpp" 75 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
77 RCP<ParameterList> validParamList = rcp(
new ParameterList());
79 validParamList->set< RCP<const FactoryBase> >(
"A", Teuchos::null,
"Generating factory of the matrix A");
80 validParamList->set< RCP<const FactoryBase> >(
"Aggregates", Teuchos::null,
"Generating factory of the aggregates");
81 validParamList->set< RCP<const FactoryBase> >(
"Nullspace", Teuchos::null,
"Generating factory of the nullspace");
82 validParamList->set< RCP<const FactoryBase> >(
"UnAmalgamationInfo", Teuchos::null,
"Generating factory of UnAmalgamationInfo");
83 validParamList->set< RCP<const FactoryBase> >(
"CoarseMap", Teuchos::null,
"Generating factory of the coarse map");
84 return validParamList;
87 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
89 Input(fineLevel,
"A");
90 Input(fineLevel,
"Aggregates");
91 Input(fineLevel,
"Nullspace");
92 Input(fineLevel,
"UnAmalgamationInfo");
93 Input(fineLevel,
"CoarseMap");
96 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
98 return BuildP(fineLevel, coarseLevel);
101 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
105 RCP<Matrix> A = Get< RCP<Matrix> > (fineLevel,
"A");
106 RCP<Aggregates> aggregates = Get< RCP<Aggregates> > (fineLevel,
"Aggregates");
107 RCP<AmalgamationInfo> amalgInfo = Get< RCP<AmalgamationInfo> >(fineLevel,
"UnAmalgamationInfo");
108 RCP<MultiVector> fineNullspace = Get< RCP<MultiVector> > (fineLevel,
"Nullspace");
109 RCP<const Map> coarseMap = Get< RCP<const Map> > (fineLevel,
"CoarseMap");
111 RCP<Matrix> Ptentative;
112 RCP<MultiVector> coarseNullspace;
113 if (!aggregates->AggregatesCrossProcessors())
114 BuildPuncoupled(A, aggregates, amalgInfo, fineNullspace, coarseMap, Ptentative, coarseNullspace);
116 BuildPcoupled (A, aggregates, amalgInfo, fineNullspace, coarseMap, Ptentative, coarseNullspace);
126 if (A->IsView(
"stridedMaps") ==
true)
127 Ptentative->CreateView(
"stridedMaps", A->getRowMap(
"stridedMaps"), coarseMap);
129 Ptentative->CreateView(
"stridedMaps", Ptentative->getRangeMap(), coarseMap);
131 Set(coarseLevel,
"Nullspace", coarseNullspace);
132 Set(coarseLevel,
"P", Ptentative);
135 RCP<ParameterList> params = rcp(
new ParameterList());
136 params->set(
"printLoadBalancingInfo",
true);
141 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
143 BuildPuncoupled(RCP<Matrix> A, RCP<Aggregates> aggregates, RCP<AmalgamationInfo> amalgInfo, RCP<MultiVector> fineNullspace,
144 RCP<const Map> coarseMap, RCP<Matrix>& Ptentative, RCP<MultiVector>& coarseNullspace)
const {
145 RCP<const Map> rowMap = A->getRowMap();
146 RCP<const Map> colMap = A->getColMap();
148 const size_t numRows = rowMap->getNodeNumElements();
150 typedef Teuchos::ScalarTraits<SC> STS;
151 typedef typename STS::magnitudeType Magnitude;
152 const SC zero = STS::zero();
153 const SC one = STS::one();
154 const LO INVALID = Teuchos::OrdinalTraits<LO>::invalid();
156 const GO numAggs = aggregates->GetNumAggregates();
157 const size_t NSDim = fineNullspace->getNumVectors();
162 bool goodMap = isGoodMap(*rowMap, *colMap);
168 ArrayRCP<LO> aggStart;
169 ArrayRCP<LO> aggToRowMapLO;
170 ArrayRCP<GO> aggToRowMapGO;
172 amalgInfo->UnamalgamateAggregatesLO(*aggregates, aggStart, aggToRowMapLO);
173 GetOStream(
Runtime1) <<
"Column map is consistent with the row map, good." << std::endl;
176 amalgInfo->UnamalgamateAggregates(*aggregates, aggStart, aggToRowMapGO);
177 GetOStream(
Warnings0) <<
"Column map is not consistent with the row map\n" 178 <<
"using GO->LO conversion with performance penalty" << std::endl;
181 coarseNullspace = MultiVectorFactory::Build(coarseMap, NSDim);
184 ArrayRCP<ArrayRCP<const SC> > fineNS (NSDim);
185 ArrayRCP<ArrayRCP<SC> > coarseNS(NSDim);
186 for (
size_t i = 0; i < NSDim; i++) {
187 fineNS[i] = fineNullspace->getData(i);
188 if (coarseMap->getNodeNumElements() > 0)
189 coarseNS[i] = coarseNullspace->getDataNonConst(i);
192 size_t nnzEstimate = numRows * NSDim;
196 RCP<CrsMatrix> PtentCrs = rcp_dynamic_cast<CrsMatrixWrap>(Ptentative)->getCrsMatrix();
198 ArrayRCP<size_t> iaPtent;
199 ArrayRCP<LO> jaPtent;
200 ArrayRCP<SC> valPtent;
202 PtentCrs->allocateAllValues(nnzEstimate, iaPtent, jaPtent, valPtent);
204 ArrayView<size_t> ia = iaPtent();
205 ArrayView<LO> ja = jaPtent();
206 ArrayView<SC> val = valPtent();
209 for (
size_t i = 1; i <= numRows; i++)
210 ia[i] = ia[i-1] + NSDim;
212 for (
size_t j = 0; j < nnzEstimate; j++) {
217 for (
GO agg = 0; agg < numAggs; agg++) {
218 LO aggSize = aggStart[agg+1] - aggStart[agg];
225 Teuchos::SerialDenseMatrix<LO,SC> localQR(aggSize, NSDim);
227 for (
size_t j = 0; j < NSDim; j++)
228 for (
LO k = 0; k < aggSize; k++)
229 localQR(k,j) = fineNS[j][aggToRowMapLO[aggStart[agg]+k]];
231 for (
size_t j = 0; j < NSDim; j++)
232 for (
LO k = 0; k < aggSize; k++)
233 localQR(k,j) = fineNS[j][rowMap->getLocalElement(aggToRowMapGO[aggStart[agg]+k])];
237 for (
size_t j = 0; j < NSDim; j++) {
238 bool bIsZeroNSColumn =
true;
240 for (
LO k = 0; k < aggSize; k++)
241 if (localQR(k,j) != zero)
242 bIsZeroNSColumn =
false;
245 "MueLu::TentativePFactory::MakeTentative: fine level NS part has a zero column");
250 if (aggSize >= Teuchos::as<LO>(NSDim)) {
254 Magnitude norm = STS::magnitude(zero);
255 for (
size_t k = 0; k < Teuchos::as<size_t>(aggSize); k++)
256 norm += STS::magnitude(localQR(k,0)*localQR(k,0));
257 norm = Teuchos::ScalarTraits<Magnitude>::squareroot(norm);
260 coarseNS[0][offset] = norm;
263 for (
LO i = 0; i < aggSize; i++)
264 localQR(i,0) /= norm;
267 Teuchos::SerialQRDenseSolver<LO,SC> qrSolver;
268 qrSolver.setMatrix(Teuchos::rcp(&localQR,
false));
272 for (
size_t j = 0; j < NSDim; j++)
273 for (
size_t k = 0; k <= j; k++)
274 coarseNS[j][offset+k] = localQR(k,j);
279 Teuchos::RCP<Teuchos::SerialDenseMatrix<LO,SC> > qFactor = qrSolver.getQ();
280 for (
size_t j = 0; j < NSDim; j++)
281 for (
size_t i = 0; i < Teuchos::as<size_t>(aggSize); i++)
282 localQR(i,j) = (*qFactor)(i,j);
313 for (
size_t j = 0; j < NSDim; j++)
314 for (
size_t k = 0; k < NSDim; k++)
315 if (k < as<size_t>(aggSize))
316 coarseNS[j][offset+k] = localQR(k,j);
318 coarseNS[j][offset+k] = (k == j ? one : zero);
321 for (
size_t i = 0; i < as<size_t>(aggSize); i++)
322 for (
size_t j = 0; j < NSDim; j++)
323 localQR(i,j) = (j == i ? one : zero);
328 for (
LO j = 0; j < aggSize; j++) {
329 LO localRow = (goodMap ? aggToRowMapLO[aggStart[agg]+j] : rowMap->getLocalElement(aggToRowMapGO[aggStart[agg]+j]));
331 size_t rowStart = ia[localRow];
332 for (
size_t k = 0, lnnz = 0; k < NSDim; k++) {
334 if (localQR(j,k) != zero) {
335 ja [rowStart+lnnz] = offset + k;
336 val[rowStart+lnnz] = localQR(j,k);
345 size_t ia_tmp = 0, nnz = 0;
346 for (
size_t i = 0; i < numRows; i++) {
347 for (
size_t j = ia_tmp; j < ia[i+1]; j++)
348 if (ja[j] != INVALID) {
360 jaPtent .resize(nnz);
361 valPtent.resize(nnz);
364 GetOStream(
Runtime1) <<
"TentativePFactory : aggregates do not cross process boundaries" << std::endl;
366 PtentCrs->setAllValues(iaPtent, jaPtent, valPtent);
367 PtentCrs->expertStaticFillComplete(coarseMap, A->getDomainMap());
370 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
372 BuildPcoupled(RCP<Matrix> A, RCP<Aggregates> aggregates, RCP<AmalgamationInfo> amalgInfo, RCP<MultiVector> fineNullspace,
373 RCP<const Map> coarseMap, RCP<Matrix>& Ptentative, RCP<MultiVector>& coarseNullspace)
const {
374 typedef Teuchos::ScalarTraits<SC> STS;
375 typedef typename STS::magnitudeType Magnitude;
376 const SC zero = STS::zero();
377 const SC one = STS::one();
380 GO numAggs = aggregates->GetNumAggregates();
386 ArrayRCP<LO> aggStart;
387 ArrayRCP< GO > aggToRowMap;
388 amalgInfo->UnamalgamateAggregates(*aggregates, aggStart, aggToRowMap);
392 for (
GO i=0; i<numAggs; ++i) {
393 LO sizeOfThisAgg = aggStart[i+1] - aggStart[i];
394 if (sizeOfThisAgg > maxAggSize) maxAggSize = sizeOfThisAgg;
398 const size_t NSDim = fineNullspace->getNumVectors();
401 GO indexBase=A->getRowMap()->getIndexBase();
403 const RCP<const Map> nonUniqueMap = amalgInfo->ComputeUnamalgamatedImportDofMap(*aggregates);
404 const RCP<const Map> uniqueMap = A->getDomainMap();
405 RCP<const Import> importer = ImportFactory::Build(uniqueMap, nonUniqueMap);
406 RCP<MultiVector> fineNullspaceWithOverlap = MultiVectorFactory::Build(nonUniqueMap,NSDim);
407 fineNullspaceWithOverlap->doImport(*fineNullspace,*importer,
Xpetra::INSERT);
410 ArrayRCP< ArrayRCP<const SC> > fineNS(NSDim);
411 for (
size_t i=0; i<NSDim; ++i)
412 fineNS[i] = fineNullspaceWithOverlap->getData(i);
415 coarseNullspace = MultiVectorFactory::Build(coarseMap, NSDim);
417 ArrayRCP< ArrayRCP<SC> > coarseNS(NSDim);
418 for (
size_t i=0; i<NSDim; ++i)
419 if (coarseMap->getNodeNumElements() > 0) coarseNS[i] = coarseNullspace->getDataNonConst(i);
424 RCP<const Map > rowMapForPtent = A->getRowMap();
425 const Map& rowMapForPtentRef = *rowMapForPtent;
429 RCP<const Map> colMap = A->getColMap();
431 RCP<const Map > ghostQMap;
432 RCP<MultiVector> ghostQvalues;
433 Array<RCP<Xpetra::Vector<GO,LO,GO,Node> > > ghostQcolumns;
434 RCP<Xpetra::Vector<GO,LO,GO,Node> > ghostQrowNums;
435 ArrayRCP< ArrayRCP<SC> > ghostQvals;
436 ArrayRCP< ArrayRCP<GO> > ghostQcols;
437 ArrayRCP< GO > ghostQrows;
440 for (
LO j=0; j<numAggs; ++j) {
441 for (
LO k=aggStart[j]; k<aggStart[j+1]; ++k) {
442 if (rowMapForPtentRef.isNodeGlobalElement(aggToRowMap[k]) ==
false) {
443 ghostGIDs.push_back(aggToRowMap[k]);
447 ghostQMap = MapFactory::Build(A->getRowMap()->lib(),
448 Teuchos::OrdinalTraits<Xpetra::global_size_t>::invalid(),
450 indexBase, A->getRowMap()->getComm());
452 ghostQvalues = MultiVectorFactory::Build(ghostQMap,NSDim);
456 ghostQcolumns.resize(NSDim);
457 for (
size_t i=0; i<NSDim; ++i)
460 if (ghostQvalues->getLocalLength() > 0) {
461 ghostQvals.resize(NSDim);
462 ghostQcols.resize(NSDim);
463 for (
size_t i=0; i<NSDim; ++i) {
464 ghostQvals[i] = ghostQvalues->getDataNonConst(i);
465 ghostQcols[i] = ghostQcolumns[i]->getDataNonConst(0);
467 ghostQrows = ghostQrowNums->getDataNonConst(0);
471 importer = ImportFactory::Build(ghostQMap, A->getRowMap());
474 Teuchos::SerialQRDenseSolver<LO,SC> qrSolver;
477 Array<GO> globalColPtr(maxAggSize*NSDim,0);
478 Array<LO> localColPtr(maxAggSize*NSDim,0);
479 Array<SC> valPtr(maxAggSize*NSDim,0.);
482 const Map& coarseMapRef = *coarseMap;
485 ArrayRCP<size_t> ptent_rowptr;
486 ArrayRCP<LO> ptent_colind;
487 ArrayRCP<Scalar> ptent_values;
490 ArrayView<size_t> rowptr_v;
491 ArrayView<LO> colind_v;
492 ArrayView<Scalar> values_v;
495 Array<size_t> rowptr_temp;
496 Array<LO> colind_temp;
497 Array<Scalar> values_temp;
499 RCP<CrsMatrix> PtentCrs;
501 RCP<CrsMatrixWrap> PtentCrsWrap = rcp(
new CrsMatrixWrap(rowMapForPtent, NSDim,
Xpetra::StaticProfile));
502 PtentCrs = PtentCrsWrap->getCrsMatrix();
503 Ptentative = PtentCrsWrap;
509 const Map& nonUniqueMapRef = *nonUniqueMap;
511 size_t total_nnz_count=0;
513 for (
GO agg=0; agg<numAggs; ++agg)
515 LO myAggSize = aggStart[agg+1]-aggStart[agg];
518 Teuchos::SerialDenseMatrix<LO,SC> localQR(myAggSize, NSDim);
519 for (
size_t j=0; j<NSDim; ++j) {
520 bool bIsZeroNSColumn =
true;
521 for (
LO k=0; k<myAggSize; ++k)
526 SC nsVal = fineNS[j][ nonUniqueMapRef.getLocalElement(aggToRowMap[aggStart[agg]+k]) ];
527 localQR(k,j) = nsVal;
528 if (nsVal != zero) bIsZeroNSColumn =
false;
531 GetOStream(
Runtime1,-1) <<
"length of fine level nsp: " << fineNullspace->getGlobalLength() << std::endl;
532 GetOStream(
Runtime1,-1) <<
"length of fine level nsp w overlap: " << fineNullspaceWithOverlap->getGlobalLength() << std::endl;
533 GetOStream(
Runtime1,-1) <<
"(local?) aggId=" << agg << std::endl;
534 GetOStream(
Runtime1,-1) <<
"aggSize=" << myAggSize << std::endl;
535 GetOStream(
Runtime1,-1) <<
"agg DOF=" << k << std::endl;
536 GetOStream(
Runtime1,-1) <<
"NS vector j=" << j << std::endl;
537 GetOStream(
Runtime1,-1) <<
"j*myAggSize + k = " << j*myAggSize + k << std::endl;
538 GetOStream(
Runtime1,-1) <<
"aggToRowMap["<<agg<<
"][" << k <<
"] = " << aggToRowMap[aggStart[agg]+k] << std::endl;
539 GetOStream(
Runtime1,-1) <<
"id aggToRowMap[agg][k]=" << aggToRowMap[aggStart[agg]+k] <<
" is global element in nonUniqueMap = " <<
540 nonUniqueMapRef.isNodeGlobalElement(aggToRowMap[aggStart[agg]+k]) << std::endl;
541 GetOStream(
Runtime1,-1) <<
"colMap local id aggToRowMap[agg][k]=" << nonUniqueMapRef.getLocalElement(aggToRowMap[aggStart[agg]+k]) << std::endl;
542 GetOStream(
Runtime1,-1) <<
"fineNS...=" << fineNS[j][ nonUniqueMapRef.getLocalElement(aggToRowMap[aggStart[agg]+k]) ] << std::endl;
543 GetOStream(
Errors,-1) <<
"caught an error!" << std::endl;
546 TEUCHOS_TEST_FOR_EXCEPTION(bIsZeroNSColumn ==
true,
Exceptions::RuntimeError,
"MueLu::TentativePFactory::MakeTentative: fine level NS part has a zero column. Error.");
551 if(myAggSize >= Teuchos::as<LocalOrdinal>(NSDim)) {
556 SC tau = localQR(0,0);
561 for (
size_t k = 0; k < Teuchos::as<size_t>(myAggSize); ++k) {
562 Magnitude tmag = STS::magnitude(localQR(k,0));
565 dtemp = Teuchos::ScalarTraits<Magnitude>::squareroot(dtemp);
567 localQR(0,0) = dtemp;
569 qrSolver.setMatrix( Teuchos::rcp(&localQR,
false) );
576 for (
size_t j=0; j<NSDim; ++j) {
577 for (
size_t k=0; k<=j; ++k) {
579 if (coarseMapRef.isNodeLocalElement(offset+k)) {
580 coarseNS[j][offset+k] = localQR(k, j);
584 GetOStream(
Errors,-1) <<
"caught error in coarseNS insert, j="<<j<<
", offset+k = "<<offset+k<<std::endl;
594 Magnitude dtemp = Teuchos::ScalarTraits<SC>::magnitude(localQR(0,0));
597 for (LocalOrdinal i=0; i<myAggSize; ++i) {
598 localQR(i,0) *= dtemp ;
602 Teuchos::RCP<Teuchos::SerialDenseMatrix<LO,SC> > qFactor = qrSolver.getQ();
603 for (
size_t j=0; j<NSDim; j++) {
604 for (
size_t i = 0; i < Teuchos::as<size_t>(myAggSize); i++) {
605 localQR(i,j) = (*qFactor)(i,j);
615 for (
size_t j = 0; j < NSDim; j++)
616 for (
size_t k = 0; k < NSDim; k++) {
618 "Caught error in coarseNS insert, j=" << j <<
", offset+k = " << offset+k);
620 if (k < as<size_t>(myAggSize))
621 coarseNS[j][offset+k] = localQR(k,j);
623 coarseNS[j][offset+k] = (k == j ? one : zero);
627 for (
size_t i = 0; i < as<size_t>(myAggSize); i++)
628 for (
size_t j = 0; j < NSDim; j++)
629 localQR(i,j) = (j == i ? one : zero);
636 for (
GO j=0; j<myAggSize; ++j) {
640 GO globalRow = aggToRowMap[aggStart[agg]+j];
643 if (rowMapForPtentRef.isNodeGlobalElement(globalRow) == false ) {
644 ghostQrows[qctr] = globalRow;
645 for (
size_t k=0; k<NSDim; ++k) {
646 ghostQcols[k][qctr] = coarseMapRef.getGlobalElement(agg*NSDim+k);
647 ghostQvals[k][qctr] = localQR(j,k);
652 for (
size_t k=0; k<NSDim; ++k) {
654 if (localQR(j,k) != Teuchos::ScalarTraits<SC>::zero()) {
655 localColPtr[nnz] = agg * NSDim + k;
656 globalColPtr[nnz] = coarseMapRef.getGlobalElement(localColPtr[nnz]);
657 valPtr[nnz] = localQR(j,k);
663 GetOStream(
Errors,-1) <<
"caught error in colPtr/valPtr insert, current index="<<nnz<<std::endl;
668 Ptentative->insertGlobalValues(globalRow,globalColPtr.view(0,nnz),valPtr.view(0,nnz));
671 GetOStream(
Errors,-1) <<
"pid " << A->getRowMap()->getComm()->getRank()
672 <<
"caught error during Ptent row insertion, global row " 673 << globalRow << std::endl;
684 GetOStream(
Runtime1) <<
"TentativePFactory : aggregates may cross process boundaries" << std::endl;
688 targetQrowNums->putScalar(-1);
689 targetQrowNums->doImport(*ghostQrowNums,*importer,
Xpetra::INSERT);
690 ArrayRCP< GO > targetQrows = targetQrowNums->getDataNonConst(0);
693 Array<GO> gidsToImport;
694 gidsToImport.reserve(targetQrows.size());
695 for (
typename ArrayRCP<GO>::iterator r=targetQrows.begin(); r!=targetQrows.end(); ++r) {
697 gidsToImport.push_back(*r);
700 RCP<const Map > reducedMap = MapFactory::Build( A->getRowMap()->lib(),
701 Teuchos::OrdinalTraits<Xpetra::global_size_t>::invalid(),
702 gidsToImport, indexBase, A->getRowMap()->getComm() );
705 importer = ImportFactory::Build(ghostQMap, reducedMap);
707 Array<RCP<Xpetra::Vector<GO,LO,GO,Node> > > targetQcolumns(NSDim);
708 for (
size_t i=0; i<NSDim; ++i) {
710 targetQcolumns[i]->doImport(*(ghostQcolumns[i]),*importer,
Xpetra::INSERT);
712 RCP<MultiVector> targetQvalues = MultiVectorFactory::Build(reducedMap,NSDim);
715 ArrayRCP< ArrayRCP<SC> > targetQvals;
716 ArrayRCP<ArrayRCP<GO> > targetQcols;
717 if (targetQvalues->getLocalLength() > 0) {
718 targetQvals.resize(NSDim);
719 targetQcols.resize(NSDim);
720 for (
size_t i=0; i<NSDim; ++i) {
721 targetQvals[i] = targetQvalues->getDataNonConst(i);
722 targetQcols[i] = targetQcolumns[i]->getDataNonConst(0);
726 valPtr = Array<SC>(NSDim,0.);
727 globalColPtr = Array<GO>(NSDim,0);
728 for (
typename Array<GO>::iterator r=gidsToImport.begin(); r!=gidsToImport.end(); ++r) {
729 if (targetQvalues->getLocalLength() > 0) {
730 for (
size_t j=0; j<NSDim; ++j) {
731 valPtr[j] = targetQvals[j][reducedMap->getLocalElement(*r)];
732 globalColPtr[j] = targetQcols[j][reducedMap->getLocalElement(*r)];
734 Ptentative->insertGlobalValues(*r, globalColPtr.view(0,NSDim), valPtr.view(0,NSDim));
738 Ptentative->fillComplete(coarseMap, A->getDomainMap());
741 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
743 ArrayView<const GO> rowElements = rowMap.getNodeElementList();
744 ArrayView<const GO> colElements = colMap.getNodeElementList();
746 const size_t numElements = rowElements.size();
749 for (
size_t i = 0; i < numElements; i++)
750 if (rowElements[i] != colElements[i]) {
762 #define MUELU_TENTATIVEPFACTORY_SHORT 763 #endif // MUELU_TENTATIVEPFACTORY_DEF_HPP Important warning messages (one line)
Timer to be used in factories. Similar to Monitor but with additional timers.
void BuildP(Level &fineLevel, Level &coarseLevel) const
Abstract Build method.
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.
Namespace for MueLu classes and methods.
void BuildPuncoupled(RCP< Matrix > A, RCP< Aggregates > aggregates, RCP< AmalgamationInfo > amalgInfo, RCP< MultiVector > fineNullspace, RCP< const Map > coarseMap, RCP< Matrix > &Ptentative, RCP< MultiVector > &coarseNullspace) const
void BuildPcoupled(RCP< Matrix > A, RCP< Aggregates > aggregates, RCP< AmalgamationInfo > amalgInfo, RCP< MultiVector > fineNullspace, RCP< const Map > coarseMap, RCP< Matrix > &Ptentative, RCP< MultiVector > &coarseNullspace) const
Class that holds all level-specific information.
static std::string PrintMatrixInfo(const Matrix &A, const std::string &msgTag, RCP< const Teuchos::ParameterList > params=Teuchos::null)
static RCP< Vector > Build(const Teuchos::RCP< const Map > &map, bool zeroOut=true)
void Build(Level &fineLevel, Level &coarseLevel) const
Build an object with this factory.
Exception throws to report errors in the internal logical of the program.
Description of what is happening (more verbose)
void DeclareInput(Level &fineLevel, Level &coarseLevel) const
Input.
bool isGoodMap(const Map &rowMap, const Map &colMap) const