46 #ifndef MUELU_SEMICOARSENPFACTORY_DEF_HPP
47 #define MUELU_SEMICOARSENPFACTORY_DEF_HPP
53 #include <Xpetra_CrsMatrixWrap.hpp>
56 #include <Xpetra_MultiVectorFactory.hpp>
67 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
71 #define SET_VALID_ENTRY(name) validParamList->setEntry(name, MasterList::getEntry(name))
76 #undef SET_VALID_ENTRY
81 validParamList->
set<
RCP<const FactoryBase> >(
"LineDetection_VertLineIds", Teuchos::null,
"Generating factory for LineDetection vertical line ids");
82 validParamList->
set<
RCP<const FactoryBase> >(
"LineDetection_Layers", Teuchos::null,
"Generating factory for LineDetection layer ids");
83 validParamList->
set<
RCP<const FactoryBase> >(
"CoarseNumZLayers", Teuchos::null,
"Generating factory for number of coarse z-layers");
85 return validParamList;
88 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
90 Input(fineLevel,
"A");
91 Input(fineLevel,
"Nullspace");
93 Input(fineLevel,
"LineDetection_VertLineIds");
94 Input(fineLevel,
"LineDetection_Layers");
95 Input(fineLevel,
"CoarseNumZLayers");
103 bTransferCoordinates_ =
true;
105 }
else if (bTransferCoordinates_ ==
true) {
110 if (myCoordsFact == Teuchos::null) {
119 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
121 return BuildP(fineLevel, coarseLevel);
124 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
130 RCP<MultiVector> fineNullspace = Get<RCP<MultiVector> >(fineLevel,
"Nullspace");
134 LO CoarsenRate = as<LO>(pL.get<
int>(
"semicoarsen: coarsen rate"));
135 bool buildRestriction = pL.get<
bool>(
"semicoarsen: calculate nonsym restriction");
136 bool doLinear = pL.get<
bool>(
"semicoarsen: piecewise linear");
139 LO BlkSize = A->GetFixedBlockSize();
141 LO Ndofs = rowMap->getLocalNumElements();
142 LO Nnodes = Ndofs / BlkSize;
145 LO FineNumZLayers = Get<LO>(fineLevel,
"CoarseNumZLayers");
146 Teuchos::ArrayRCP<LO> TVertLineId = Get<Teuchos::ArrayRCP<LO> >(fineLevel,
"LineDetection_VertLineIds");
149 LO* LayerId = TLayerId.getRawPtr();
155 GO Ncoarse = MakeSemiCoarsenP(Nnodes, FineNumZLayers, CoarsenRate, LayerId, VertLineId,
156 BlkSize, A, P, theCoarseMap, fineNullspace, coarseNullspace, R, buildRestriction, doLinear);
159 if (A->IsView(
"stridedMaps") ==
true)
160 P->CreateView(
"stridedMaps", A->getRowMap(
"stridedMaps"), theCoarseMap);
162 P->CreateView(
"stridedMaps", P->getRangeMap(), theCoarseMap);
164 if (buildRestriction) {
165 if (A->IsView(
"stridedMaps") ==
true)
166 R->CreateView(
"stridedMaps", theCoarseMap, A->getRowMap(
"stridedMaps"));
168 R->CreateView(
"stridedMaps", theCoarseMap, R->getDomainMap());
170 if (pL.get<
bool>(
"semicoarsen: piecewise constant")) {
172 RevertToPieceWiseConstant(P, BlkSize);
174 if (pL.get<
bool>(
"semicoarsen: piecewise linear")) {
183 LO CoarseNumZLayers = FineNumZLayers * Ncoarse;
184 CoarseNumZLayers /= Ndofs;
188 Set(coarseLevel,
"P", P);
189 if (buildRestriction) Set(coarseLevel,
"RfromPfactory", R);
191 Set(coarseLevel,
"Nullspace", coarseNullspace);
194 if (bTransferCoordinates_) {
197 if (fineLevel.GetLevelID() == 0 &&
202 if (myCoordsFact == Teuchos::null) {
203 myCoordsFact = fineLevel.GetFactoryManager()->GetFactory(
"Coordinates");
205 if (fineLevel.IsAvailable(
"Coordinates", myCoordsFact.
get())) {
206 fineCoords = fineLevel.Get<
RCP<xdMV> >(
"Coordinates", myCoordsFact.
get());
220 for (
auto it = z.
begin(); it != z.
end(); ++it) {
221 if (*it > zval_max) zval_max = *it;
222 if (*it < zval_min) zval_min = *it;
225 LO myCoarseZLayers = Teuchos::as<LO>(CoarseNumZLayers);
228 if (myCoarseZLayers == 1) {
229 myZLayerCoords[0] = zval_min;
232 for (LO k = 0; k < myCoarseZLayers; ++k) {
233 myZLayerCoords[k] = k * dz;
241 LO numVertLines = Nnodes / FineNumZLayers;
242 LO numLocalCoarseNodes = numVertLines * myCoarseZLayers;
245 MapFactory::Build(fineCoords->getMap()->lib(),
247 Teuchos::as<size_t>(numLocalCoarseNodes),
248 fineCoords->getMap()->getIndexBase(),
249 fineCoords->getMap()->getComm());
251 coarseCoords->putScalar(-1.0);
257 LO cntCoarseNodes = 0;
258 for (LO vt = 0; vt < TVertLineId.size(); ++vt) {
260 LO curVertLineId = TVertLineId[vt];
262 if (cx[curVertLineId * myCoarseZLayers] == -1.0 &&
263 cy[curVertLineId * myCoarseZLayers] == -1.0) {
265 for (LO n = 0; n < myCoarseZLayers; ++n) {
266 cx[curVertLineId * myCoarseZLayers + n] = x[vt];
267 cy[curVertLineId * myCoarseZLayers + n] = y[vt];
268 cz[curVertLineId * myCoarseZLayers + n] = myZLayerCoords[n];
270 cntCoarseNodes += myCoarseZLayers;
274 if (cntCoarseNodes == numLocalCoarseNodes)
break;
278 Set(coarseLevel,
"Coordinates", coarseCoords);
282 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
313 NCpts = (
LO)ceil(temp);
315 NCpts = (
LO)floor(temp);
319 if (NCpts < 1) NCpts = 1;
329 for (i = 1; i <= NCpts; i++) {
330 (*LayerCpts)[i] = (
LO)floor(di);
337 #define MaxHorNeighborNodes 75
339 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
395 int NLayers, NVertLines, MaxNnz, NCLayers, MyLine, MyLayer;
396 int *InvLineLayer = NULL, *CptLayers = NULL, StartLayer, NStencilNodes;
397 int BlkRow = -1, dof_j, node_k, *Sub2FullMap = NULL, RowLeng;
398 int i, j, iii, col, count, index, loc, PtRow, PtCol;
399 SC *BandSol = NULL, *BandMat = NULL, TheSum, *RestrictBandMat = NULL, *RestrictBandSol = NULL;
400 int *IPIV = NULL, KL, KU, KLU, N, NRHS, LDAB, INFO;
407 int MaxStencilSize, MaxNnzPerRow;
409 int CurRow, LastGuy = -1, NewPtr, RLastGuy = -1;
412 LO *Layerdofs = NULL, *Col2Dof = NULL;
427 Nghost = Amat->getColMap()->getLocalNumElements() - Amat->getDomainMap()->getLocalNumElements();
428 if (Nghost < 0) Nghost = 0;
438 for (i = 0; i < Ntotal * DofsPerNode; i++)
439 valptr[i] = LayerId[i / DofsPerNode];
443 importer = Amat->getCrsGraph()->getImporter();
444 if (importer == Teuchos::null) {
445 importer = ImportFactory::Build(Amat->getDomainMap(), Amat->getColMap());
449 valptr = dtemp->getDataNonConst(0);
450 for (i = 0; i < Ntotal * DofsPerNode + Nghost; i++) Layerdofs[i] = valptr[i];
451 valptr = localdtemp->getDataNonConst(0);
452 for (i = 0; i < Ntotal * DofsPerNode; i++) valptr[i] = i % DofsPerNode;
456 valptr = dtemp->getDataNonConst(0);
457 for (i = 0; i < Ntotal * DofsPerNode + Nghost; i++) Col2Dof[i] = valptr[i];
461 NLayers = LayerId[0];
462 NVertLines = VertLineId[0];
468 for (i = 1; i < Ntotal; i++) {
469 if (VertLineId[i] > NVertLines) NVertLines = VertLineId[i];
470 if (LayerId[i] > NLayers) NLayers = LayerId[i];
481 InvLineLayer = TInvLineLayer.
getRawPtr();
482 for (i = 0; i < Ntotal; i++) {
483 InvLineLayer[VertLineId[i] + 1 + LayerId[i] * NVertLines] = i;
492 NCLayers = FindCpts(nz, CoarsenRate, 0, &CptLayers);
504 MaxStencilSize = CptLayers[2];
506 for (i = 3; i <= NCLayers; i++) {
507 if (MaxStencilSize < CptLayers[i] - CptLayers[i - 2])
508 MaxStencilSize = CptLayers[i] - CptLayers[i - 2];
511 if (MaxStencilSize < nz - CptLayers[NCLayers - 1] + 1)
512 MaxStencilSize = nz - CptLayers[NCLayers - 1] + 1;
524 Teuchos::ArrayRCP<SC> TBandSol = Teuchos::arcp<SC>((MaxStencilSize + 1) * DofsPerNode * DofsPerNode);
527 if (buildRestriction) {
528 TResBandSol = Teuchos::arcp<SC>((MaxStencilSize + 1) * DofsPerNode * DofsPerNode);
529 RestrictBandSol = TResBandSol.
getRawPtr();
534 KL = 2 * DofsPerNode - 1;
535 KU = 2 * DofsPerNode - 1;
537 LDAB = 2 * KL + KU + 1;
545 if (buildRestriction) {
546 TResBandMat = Teuchos::arcp<SC>(LDAB * MaxStencilSize * DofsPerNode + 1);
547 RestrictBandMat = TResBandMat.
getRawPtr();
557 Ndofs = DofsPerNode * Ntotal;
558 MaxNnz = 2 * DofsPerNode * Ndofs;
563 std::vector<size_t> stridingInfo_;
564 stridingInfo_.push_back(DofsPerNode);
567 coarseMap = StridedMapFactory::Build(rowMap->lib(),
569 NCLayers * NVertLines * DofsPerNode,
578 P =
rcp(
new CrsMatrixWrap(rowMap, coarseMap, 0));
579 coarseNullspace = MultiVectorFactory::Build(coarseMap, fineNullspace->getNumVectors());
595 LO RmaxNnz = -1, RmaxPerRow = -1;
596 if (buildRestriction) {
597 RmaxPerRow = ((
LO)ceil(((
double)(MaxNnz + 7)) / ((double)(coarseMap->getLocalNumElements()))));
598 RmaxNnz = RmaxPerRow * coarseMap->getLocalNumElements();
599 TRvals = Teuchos::arcp<SC>(1 + RmaxNnz);
601 TRptr = Teuchos::arcp<size_t>((2 + coarseMap->getLocalNumElements()));
603 TRcols = Teuchos::arcp<LO>(1 + RmaxNnz);
615 for (i = 1; i <= MaxNnz; i++) Pcols[i] = -1;
617 for (i = 1; i <= DofsPerNode * Ntotal + 1; i++) {
619 count += (2 * DofsPerNode);
621 if (buildRestriction) {
622 for (i = 1; i <= RmaxNnz; i++) Rcols[i] = -1;
624 for (i = 1; i <= (int)(coarseMap->getLocalNumElements() + 1); i++) {
638 for (MyLine = 1; MyLine <= NVertLines; MyLine += 1) {
639 for (iii = 1; iii <= NCLayers; iii += 1) {
641 MyLayer = CptLayers[iii];
655 StartLayer = CptLayers[iii - 1] + 1;
660 NStencilNodes = CptLayers[iii + 1] - StartLayer;
662 NStencilNodes = NLayers - StartLayer + 1;
664 N = NStencilNodes * DofsPerNode;
671 for (i = 0; i < NStencilNodes * DofsPerNode * DofsPerNode; i++) BandSol[i] = 0.0;
672 for (i = 0; i < LDAB * N; i++) BandMat[i] = 0.0;
673 if (buildRestriction) {
674 for (i = 0; i < NStencilNodes * DofsPerNode * DofsPerNode; i++) RestrictBandSol[i] = 0.0;
675 for (i = 0; i < LDAB * N; i++) RestrictBandMat[i] = 0.0;
684 for (node_k = 1; node_k <= NStencilNodes; node_k++) {
689 BlkRow = InvLineLayer[MyLine + (StartLayer + node_k - 2) * NVertLines] + 1;
690 Sub2FullMap[node_k] = BlkRow;
703 for (
int dof_i = 0; dof_i < DofsPerNode; dof_i++) {
704 j = (BlkRow - 1) * DofsPerNode + dof_i;
707 Amat->getLocalRowView(j, AAcols, AAvals);
710 RowLeng = AAvals.
size();
714 for (i = 0; i < RowLeng; i++) {
715 LayDiff[i] = Layerdofs[Acols[i]] - StartLayer - node_k + 2;
723 PtRow = (node_k - 1) * DofsPerNode + dof_i + 1;
724 for (dof_j = 0; dof_j < DofsPerNode; dof_j++) {
725 PtCol = (node_k - 1) * DofsPerNode + dof_j + 1;
729 for (i = 0; i < RowLeng; i++) {
730 if ((LayDiff[i] == 0) && (Col2Dof[Acols[i]] == dof_j))
733 index = LDAB * (PtCol - 1) + KLU + PtRow - PtCol;
734 BandMat[index] = TheSum;
735 if (buildRestriction) RestrictBandMat[index] = TheSum;
736 if (node_k != NStencilNodes) {
740 for (i = 0; i < RowLeng; i++) {
741 if ((LayDiff[i] == 1) && (Col2Dof[Acols[i]] == dof_j))
744 j = PtCol + DofsPerNode;
745 index = LDAB * (j - 1) + KLU + PtRow - j;
746 BandMat[index] = TheSum;
747 if (buildRestriction) RestrictBandMat[index] = TheSum;
753 for (i = 0; i < RowLeng; i++) {
754 if ((LayDiff[i] == -1) && (Col2Dof[Acols[i]] == dof_j))
757 j = PtCol - DofsPerNode;
758 index = LDAB * (j - 1) + KLU + PtRow - j;
759 BandMat[index] = TheSum;
760 if (buildRestriction) RestrictBandMat[index] = TheSum;
765 node_k = MyLayer - StartLayer + 1;
766 for (
int dof_i = 0; dof_i < DofsPerNode; dof_i++) {
769 PtRow = (node_k - 1) * DofsPerNode + dof_i + 1;
770 BandSol[(dof_i)*DofsPerNode * NStencilNodes + PtRow - 1] = 1.;
771 if (buildRestriction) RestrictBandSol[(dof_i)*DofsPerNode * NStencilNodes + PtRow - 1] = 1.;
773 for (dof_j = 0; dof_j < DofsPerNode; dof_j++) {
774 PtCol = (node_k - 1) * DofsPerNode + dof_j + 1;
775 index = LDAB * (PtCol - 1) + KLU + PtRow - PtCol;
777 BandMat[index] = 1.0;
779 BandMat[index] = 0.0;
780 if (buildRestriction) {
781 index = LDAB * (PtRow - 1) + KLU + PtCol - PtRow;
783 RestrictBandMat[index] = 1.0;
785 RestrictBandMat[index] = 0.0;
787 if (node_k != NStencilNodes) {
788 PtCol = (node_k)*DofsPerNode + dof_j + 1;
789 index = LDAB * (PtCol - 1) + KLU + PtRow - PtCol;
790 BandMat[index] = 0.0;
791 if (buildRestriction) {
792 index = LDAB * (PtRow - 1) + KLU + PtCol - PtRow;
793 RestrictBandMat[index] = 0.0;
797 PtCol = (node_k - 2) * DofsPerNode + dof_j + 1;
798 index = LDAB * (PtCol - 1) + KLU + PtRow - PtCol;
799 BandMat[index] = 0.0;
800 if (buildRestriction) {
801 index = LDAB * (PtRow - 1) + KLU + PtCol - PtRow;
802 RestrictBandMat[index] = 0.0;
810 lapack.
GBTRF(N, N, KL, KU, BandMat, LDAB, IPIV, &INFO);
814 lapack.
GBTRS(notrans[0], N, KL, KU, NRHS, BandMat, LDAB, IPIV,
819 for (dof_j = 0; dof_j < DofsPerNode; dof_j++) {
820 for (
int dof_i = 0; dof_i < DofsPerNode; dof_i++) {
821 for (i = 1; i <= NStencilNodes; i++) {
822 index = (Sub2FullMap[i] - 1) * DofsPerNode + dof_i + 1;
824 Pcols[loc] = (col - 1) * DofsPerNode + dof_j + 1;
825 Pvals[loc] = BandSol[dof_j * DofsPerNode * NStencilNodes +
826 (i - 1) * DofsPerNode + dof_i];
827 Pptr[index] = Pptr[index] + 1;
831 if (buildRestriction) {
832 lapack.
GBTRF(N, N, KL, KU, RestrictBandMat, LDAB, IPIV, &INFO);
834 lapack.
GBTRS(trans[0], N, KL, KU, NRHS, RestrictBandMat, LDAB, IPIV, RestrictBandSol, N, &INFO);
836 for (dof_j = 0; dof_j < DofsPerNode; dof_j++) {
837 for (
int dof_i = 0; dof_i < DofsPerNode; dof_i++) {
838 for (i = 1; i <= NStencilNodes; i++) {
839 index = (col - 1) * DofsPerNode + dof_j + 1;
841 Rcols[loc] = (Sub2FullMap[i] - 1) * DofsPerNode + dof_i + 1;
842 Rvals[loc] = RestrictBandSol[dof_j * DofsPerNode * NStencilNodes +
843 (i - 1) * DofsPerNode + dof_i];
844 Rptr[index] = Rptr[index] + 1;
850 int denom1 = MyLayer - StartLayer + 1;
851 int denom2 = StartLayer + NStencilNodes - MyLayer;
852 for (
int dof_i = 0; dof_i < DofsPerNode; dof_i++) {
853 for (i = 1; i <= NStencilNodes; i++) {
854 index = (InvLineLayer[MyLine + (StartLayer + i - 2) * NVertLines]) * DofsPerNode + dof_i + 1;
856 Pcols[loc] = (col - 1) * DofsPerNode + dof_i + 1;
858 Pvals[loc] = 1.0 + ((double)(denom1 - i)) / ((
double)denom2);
860 Pvals[loc] = ((double)(i)) / ((
double)denom1);
861 Pptr[index] = Pptr[index] + 1;
869 BlkRow = InvLineLayer[MyLine + (MyLayer - 1) * NVertLines] + 1;
870 for (
int k = 0; k < static_cast<int>(fineNullspace->getNumVectors()); k++) {
873 for (
int dof_i = 0; dof_i < DofsPerNode; dof_i++) {
874 OneCNull[(col - 1) * DofsPerNode + dof_i] = OneFNull[(BlkRow - 1) * DofsPerNode + dof_i];
889 for (
size_t ii = 1; ii <= Pptr[Ntotal * DofsPerNode] - 1; ii++) {
890 if (ii == Pptr[CurRow]) {
891 Pptr[CurRow] = LastGuy;
893 while (ii > Pptr[CurRow]) {
894 Pptr[CurRow] = LastGuy;
898 if (Pcols[ii] != -1) {
899 Pcols[NewPtr - 1] = Pcols[ii] - 1;
900 Pvals[NewPtr - 1] = Pvals[ii];
905 for (i = CurRow; i <= Ntotal * DofsPerNode; i++) Pptr[CurRow] = LastGuy;
910 for (i = -Ntotal * DofsPerNode + 1; i >= 2; i--) {
911 Pptr[i - 1] = Pptr[i - 2];
916 if (buildRestriction) {
919 for (
size_t ii = 1; ii <= Rptr[coarseMap->getLocalNumElements()] - 1; ii++) {
920 if (ii == Rptr[CurRow]) {
921 Rptr[CurRow] = RLastGuy;
923 while (ii > Rptr[CurRow]) {
924 Rptr[CurRow] = RLastGuy;
928 if (Rcols[ii] != -1) {
929 Rcols[NewPtr - 1] = Rcols[ii] - 1;
930 Rvals[NewPtr - 1] = Rvals[ii];
935 for (i = CurRow; i <= ((int)coarseMap->getLocalNumElements()); i++) Rptr[CurRow] = RLastGuy;
936 for (i = -coarseMap->getLocalNumElements() + 1; i >= 2; i--) {
937 Rptr[i - 1] = Rptr[i - 2];
945 RCP<CrsMatrix> PCrs = rcp_dynamic_cast<CrsMatrixWrap>(P)->getCrsMatrix();
946 LO nnz = Pptr[Ndofs];
947 PCrs->allocateAllValues(nnz, rcpRowPtr, rcpColumns, rcpValues);
955 for (
LO ii = 0; ii <= Ndofs; ii++) rowPtr[ii] = Pptr[ii];
956 for (
LO ii = 0; ii < nnz; ii++) columns[ii] = Pcols[ii];
957 for (
LO ii = 0; ii < nnz; ii++) values[ii] = Pvals[ii];
958 PCrs->setAllValues(rcpRowPtr, rcpColumns, rcpValues);
959 PCrs->expertStaticFillComplete(coarseMap, Amat->getDomainMap());
964 if (buildRestriction) {
965 R =
rcp(
new CrsMatrixWrap(coarseMap, rowMap, 0));
966 RCrs = rcp_dynamic_cast<CrsMatrixWrap>(R)->getCrsMatrix();
967 nnz = Rptr[coarseMap->getLocalNumElements()];
968 RCrs->allocateAllValues(nnz, RrcpRowPtr, RrcpColumns, RrcpValues);
976 for (
LO ii = 0; ii <= ((
LO)coarseMap->getLocalNumElements()); ii++) RrowPtr[ii] = Rptr[ii];
977 for (
LO ii = 0; ii < nnz; ii++) Rcolumns[ii] = Rcols[ii];
978 for (
LO ii = 0; ii < nnz; ii++) Rvalues[ii] = Rvals[ii];
979 RCrs->setAllValues(RrcpRowPtr, RrcpColumns, RrcpValues);
980 RCrs->expertStaticFillComplete(Amat->getRangeMap(), coarseMap);
983 return NCLayers * NVertLines * DofsPerNode;
985 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
999 for (
size_t i = 0; i < P->getRowMap()->getLocalNumElements(); i++) {
1000 P->getLocalRowView(i, inds, vals1);
1002 size_t nnz = inds.
size();
1005 LO largestIndex = -1;
1006 Scalar largestValue = ZERO;
1009 LO rowDof = i % BlkSize;
1010 for (
size_t j = 0; j < nnz; j++) {
1012 if (inds[j] % BlkSize == rowDof) {
1013 largestValue = vals[j];
1014 largestIndex = (int)j;
1019 if (largestIndex != -1)
1020 vals[largestIndex] = ONE;
1030 #define MUELU_SEMICOARSENPFACTORY_SHORT
1031 #endif // MUELU_SEMICOARSENPFACTORY_DEF_HPP
MueLu::DefaultLocalOrdinal LocalOrdinal
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.
Timer to be used in factories. Similar to Monitor but with additional timers.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
void DeclareInput(Level &fineLevel, Level &coarseLevel) const
Input.
void GBTRF(const OrdinalType &m, const OrdinalType &n, const OrdinalType &kl, const OrdinalType &ku, ScalarType *A, const OrdinalType &lda, OrdinalType *IPIV, OrdinalType *info) const
LO MakeSemiCoarsenP(LO const Ntotal, LO const nz, LO const CoarsenRate, LO const LayerId[], LO const VertLineId[], LO const DofsPerNode, RCP< Matrix > &Amat, RCP< Matrix > &P, RCP< const Map > &coarseMap, const RCP< MultiVector > fineNullspace, RCP< MultiVector > &coarseNullspace, RCP< Matrix > &R, bool buildRestriction, bool doLinear) const
void Build(Level &fineLevel, Level &coarseLevel) const
Build an object with this factory.
static const NoFactory * get()
void RevertToPieceWiseConstant(RCP< Matrix > &P, LO BlkSize) const
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
MueLu::DefaultScalar Scalar
void BuildP(Level &fineLevel, Level &coarseLevel) const
Abstract Build method.
Class that holds all level-specific information.
LO FindCpts(LO const PtsPerLine, LO const CoarsenRate, LO const Thin, LO **LayerCpts) const
void Set(const std::string &ename, const T &entry, const FactoryBase *factory=NoFactory::get())
static RCP< Vector > Build(const Teuchos::RCP< const Map > &map, bool zeroOut=true)
void GBTRS(const char &TRANS, const OrdinalType &n, const OrdinalType &kl, const OrdinalType &ku, const OrdinalType &nrhs, const ScalarType *A, const OrdinalType &lda, const OrdinalType *IPIV, ScalarType *B, const OrdinalType &ldb, OrdinalType *info) const
#define SET_VALID_ENTRY(name)
#define MaxHorNeighborNodes
const RCP< const FactoryManagerBase > GetFactoryManager()
returns the current factory manager
int GetLevelID() const
Return level number.
Exception throws to report errors in the internal logical of the program.
void DeclareInput(const std::string &ename, const FactoryBase *factory, const FactoryBase *requestedBy=NoFactory::get())
Callback from FactoryBase::CallDeclareInput() and FactoryBase::DeclareInput()
bool IsAvailable(const std::string &ename, const FactoryBase *factory=NoFactory::get()) const
Test whether a need's value has been saved.