46 #ifndef MUELU_SEMICOARSENPFACTORY_DEF_HPP 47 #define MUELU_SEMICOARSENPFACTORY_DEF_HPP 51 #include <Teuchos_LAPACK.hpp> 53 #include <Xpetra_CrsMatrixWrap.hpp> 54 #include <Xpetra_ImportFactory.hpp> 55 #include <Xpetra_Matrix.hpp> 56 #include <Xpetra_MultiVectorFactory.hpp> 57 #include <Xpetra_VectorFactory.hpp> 67 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
69 RCP<ParameterList> validParamList = rcp(
new ParameterList());
71 #define SET_VALID_ENTRY(name) validParamList->setEntry(name, MasterList::getEntry(name)) 73 #undef SET_VALID_ENTRY 74 validParamList->set< RCP<const FactoryBase> >(
"A", Teuchos::null,
"Generating factory of the matrix A");
75 validParamList->set< RCP<const FactoryBase> >(
"Nullspace", Teuchos::null,
"Generating factory of the nullspace");
76 validParamList->set< RCP<const FactoryBase> >(
"Coordinates", Teuchos::null,
"Generating factory for coorindates");
78 validParamList->set< RCP<const FactoryBase> >(
"LineDetection_VertLineIds", Teuchos::null,
"Generating factory for LineDetection information");
79 validParamList->set< RCP<const FactoryBase> >(
"LineDetection_Layers", Teuchos::null,
"Generating factory for LineDetection information");
80 validParamList->set< RCP<const FactoryBase> >(
"CoarseNumZLayers", Teuchos::null,
"Generating factory for LineDetection information");
82 return validParamList;
85 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
87 Input(fineLevel,
"A");
88 Input(fineLevel,
"Nullspace");
90 Input(fineLevel,
"LineDetection_VertLineIds");
91 Input(fineLevel,
"LineDetection_Layers");
92 Input(fineLevel,
"CoarseNumZLayers");
100 bTransferCoordinates_ =
true;
102 }
else if (bTransferCoordinates_ ==
true){
106 RCP<const FactoryBase> myCoordsFact = GetFactory(
"Coordinates");
107 if (myCoordsFact == Teuchos::null) { myCoordsFact = fineLevel.
GetFactoryManager()->GetFactory(
"Coordinates"); }
108 if (fineLevel.
IsAvailable(
"Coordinates", myCoordsFact.get())) {
109 fineLevel.
DeclareInput(
"Coordinates", myCoordsFact.get(),
this);
110 bTransferCoordinates_ =
true;
115 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
117 return BuildP(fineLevel, coarseLevel);
120 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
125 RCP<Matrix> A = Get< RCP<Matrix> > (fineLevel,
"A");
126 RCP<MultiVector> fineNullspace = Get< RCP<MultiVector> > (fineLevel,
"Nullspace");
129 const ParameterList& pL = GetParameterList();
130 LO CoarsenRate = as<LO>(pL.get<
int>(
"semicoarsen: coarsen rate"));
133 LO BlkSize = A->GetFixedBlockSize();
134 RCP<const Map> rowMap = A->getRowMap();
135 LO Ndofs = rowMap->getNodeNumElements();
136 LO Nnodes = Ndofs/BlkSize;
139 LO FineNumZLayers = Get< LO >(fineLevel,
"CoarseNumZLayers");
140 Teuchos::ArrayRCP<LO> TVertLineId = Get< Teuchos::ArrayRCP<LO> > (fineLevel,
"LineDetection_VertLineIds");
141 Teuchos::ArrayRCP<LO> TLayerId = Get< Teuchos::ArrayRCP<LO> > (fineLevel,
"LineDetection_Layers");
142 LO* VertLineId = TVertLineId.getRawPtr();
143 LO* LayerId = TLayerId.getRawPtr();
146 RCP<const Map> theCoarseMap;
148 RCP<MultiVector> coarseNullspace;
149 GO Ncoarse = MakeSemiCoarsenP(Nnodes,FineNumZLayers,CoarsenRate,LayerId,VertLineId,
150 BlkSize, A, P, theCoarseMap, fineNullspace,coarseNullspace);
153 if (A->IsView(
"stridedMaps") ==
true)
154 P->CreateView(
"stridedMaps", A->getRowMap(
"stridedMaps"), theCoarseMap);
156 P->CreateView(
"stridedMaps", P->getRangeMap(), theCoarseMap);
162 LO CoarseNumZLayers = FineNumZLayers*Ncoarse;
163 CoarseNumZLayers /= Ndofs;
167 Set(coarseLevel,
"P", P);
169 Set(coarseLevel,
"Nullspace", coarseNullspace);
172 if(bTransferCoordinates_) {
174 typedef Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::coordinateType,LO,GO,NO> xdMV;
175 RCP<xdMV> fineCoords = Teuchos::null;
176 if (fineLevel.GetLevelID() == 0 &&
178 fineCoords = fineLevel.Get< RCP<xdMV> >(
"Coordinates",
NoFactory::get());
180 RCP<const FactoryBase> myCoordsFact = GetFactory(
"Coordinates");
181 if (myCoordsFact == Teuchos::null) { myCoordsFact = fineLevel.GetFactoryManager()->GetFactory(
"Coordinates"); }
182 if (fineLevel.IsAvailable(
"Coordinates", myCoordsFact.get())) {
183 fineCoords = fineLevel.Get< RCP<xdMV> >(
"Coordinates", myCoordsFact.get());
187 TEUCHOS_TEST_FOR_EXCEPTION(fineCoords==Teuchos::null,
Exceptions::RuntimeError,
"No Coordinates found provided by the user.");
189 TEUCHOS_TEST_FOR_EXCEPTION(fineCoords->getNumVectors() != 3,
Exceptions::RuntimeError,
"Three coordinates arrays must be supplied if line detection orientation not given.");
190 ArrayRCP<typename Teuchos::ScalarTraits<Scalar>::coordinateType> x = fineCoords->getDataNonConst(0);
191 ArrayRCP<typename Teuchos::ScalarTraits<Scalar>::coordinateType> y = fineCoords->getDataNonConst(1);
192 ArrayRCP<typename Teuchos::ScalarTraits<Scalar>::coordinateType> z = fineCoords->getDataNonConst(2);
195 typename Teuchos::ScalarTraits<Scalar>::coordinateType zval_max = -Teuchos::ScalarTraits<typename Teuchos::ScalarTraits<Scalar>::coordinateType>::one() / Teuchos::ScalarTraits<typename Teuchos::ScalarTraits<Scalar>::coordinateType>::sfmin();
196 typename Teuchos::ScalarTraits<Scalar>::coordinateType zval_min = Teuchos::ScalarTraits<typename Teuchos::ScalarTraits<Scalar>::coordinateType>::one() / Teuchos::ScalarTraits<typename Teuchos::ScalarTraits<Scalar>::coordinateType>::sfmin();
197 for (
auto it = z.begin(); it != z.end(); ++it) {
198 if(*it > zval_max) zval_max = *it;
199 if(*it < zval_min) zval_min = *it;
202 LO myCoarseZLayers = Teuchos::as<LO>(CoarseNumZLayers);
204 ArrayRCP<typename Teuchos::ScalarTraits<Scalar>::coordinateType> myZLayerCoords = Teuchos::arcp<typename Teuchos::ScalarTraits<Scalar>::coordinateType>(myCoarseZLayers);
205 if(myCoarseZLayers == 1) {
206 myZLayerCoords[0] = zval_min;
208 typename Teuchos::ScalarTraits<Scalar>::coordinateType dz = (zval_max-zval_min)/(myCoarseZLayers-1);
209 for(LO k = 0; k<myCoarseZLayers; ++k) {
210 myZLayerCoords[k] = k*dz;
218 LO numVertLines = Nnodes / FineNumZLayers;
219 LO numLocalCoarseNodes = numVertLines * myCoarseZLayers;
227 RCP<const Map> coarseCoordMap =
228 MapFactory::Build (fineCoords->getMap()->lib(),
229 Teuchos::OrdinalTraits<Xpetra::global_size_t>::invalid(),
230 Teuchos::as<size_t>(numLocalCoarseNodes),
231 fineCoords->getMap()->getIndexBase(),
232 fineCoords->getMap()->getComm());
233 RCP<xdMV> coarseCoords = Xpetra::MultiVectorFactory<typename Teuchos::ScalarTraits<Scalar>::coordinateType,LO,GO,NO>::Build(coarseCoordMap, fineCoords->getNumVectors());
234 coarseCoords->putScalar(-1.0);
235 ArrayRCP<typename Teuchos::ScalarTraits<Scalar>::coordinateType> cx = coarseCoords->getDataNonConst(0);
236 ArrayRCP<typename Teuchos::ScalarTraits<Scalar>::coordinateType> cy = coarseCoords->getDataNonConst(1);
237 ArrayRCP<typename Teuchos::ScalarTraits<Scalar>::coordinateType> cz = coarseCoords->getDataNonConst(2);
240 LO cntCoarseNodes = 0;
241 for( LO vt = 0; vt < TVertLineId.size(); ++vt) {
243 LO curVertLineId = TVertLineId[vt];
245 if(cx[curVertLineId * myCoarseZLayers] == -1.0 &&
246 cy[curVertLineId * myCoarseZLayers] == -1.0) {
248 for (LO n=0; n<myCoarseZLayers; ++n) {
249 cx[curVertLineId * myCoarseZLayers + n] = x[vt];
250 cy[curVertLineId * myCoarseZLayers + n] = y[vt];
251 cz[curVertLineId * myCoarseZLayers + n] = myZLayerCoords[n];
253 cntCoarseNodes += myCoarseZLayers;
256 TEUCHOS_TEST_FOR_EXCEPTION(cntCoarseNodes > numLocalCoarseNodes,
Exceptions::RuntimeError,
"number of coarse nodes is inconsistent.");
257 if(cntCoarseNodes == numLocalCoarseNodes)
break;
261 Set(coarseLevel,
"Coordinates", coarseCoords);
266 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
290 typename Teuchos::ScalarTraits<Scalar>::coordinateType temp, RestStride, di;
295 temp = ((
typename Teuchos::ScalarTraits<Scalar>::coordinateType) (PtsPerLine+1))/((
typename Teuchos::ScalarTraits<Scalar>::coordinateType) (CoarsenRate)) - 1.0;
296 if (Thin == 1) NCpts = (LO) ceil(temp);
297 else NCpts = (LO) floor(temp);
299 TEUCHOS_TEST_FOR_EXCEPTION(PtsPerLine == 1,
Exceptions::RuntimeError,
"SemiCoarsenPFactory::FindCpts: cannot coarsen further.");
301 if (NCpts < 1) NCpts = 1;
303 FirstStride= (LO) ceil( ((
typename Teuchos::ScalarTraits<Scalar>::coordinateType) PtsPerLine+1)/( (
typename Teuchos::ScalarTraits<Scalar>::coordinateType) (NCpts+1)));
304 RestStride = ((
typename Teuchos::ScalarTraits<Scalar>::coordinateType) (PtsPerLine-FirstStride+1))/((
typename Teuchos::ScalarTraits<Scalar>::coordinateType) NCpts);
306 NCLayers = (LO) floor((((
typename Teuchos::ScalarTraits<Scalar>::coordinateType) (PtsPerLine-FirstStride+1))/RestStride)+.00001);
310 di = (
typename Teuchos::ScalarTraits<Scalar>::coordinateType) FirstStride;
311 for (i = 1; i <= NCpts; i++) {
312 (*LayerCpts)[i] = (LO) floor(di);
319 #define MaxHorNeighborNodes 75 321 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
324 LO
const VertLineId[], LO
const DofsPerNode, RCP<Matrix> & Amat, RCP<Matrix>& P,
325 RCP<const Map>& coarseMap,
const RCP<MultiVector> fineNullspace,
326 RCP<MultiVector>& coarseNullspace )
const {
378 int NLayers, NVertLines, MaxNnz, NCLayers, MyLine, MyLayer;
379 int *InvLineLayer=NULL, *CptLayers=NULL, StartLayer, NStencilNodes;
380 int BlkRow, dof_j, node_k, *Sub2FullMap=NULL, RowLeng;
381 int i, j, iii, col, count, index, loc, PtRow, PtCol;
382 SC *BandSol=NULL, *BandMat=NULL, TheSum;
383 int *IPIV=NULL, KL, KU, KLU, N, NRHS, LDAB,INFO;
387 int MaxStencilSize, MaxNnzPerRow;
389 int CurRow, LastGuy = -1, NewPtr;
392 LO *Layerdofs = NULL, *Col2Dof = NULL;
394 Teuchos::LAPACK<LO,SC> lapack;
402 Teuchos::ArrayRCP<LO> TLayDiff = Teuchos::arcp<LO>(1+MaxNnzPerRow); LayDiff = TLayDiff.getRawPtr();
404 Nghost = Amat->getColMap()->getNodeNumElements() - Amat->getDomainMap()->getNodeNumElements();
405 if (Nghost < 0) Nghost = 0;
406 Teuchos::ArrayRCP<LO> TLayerdofs= Teuchos::arcp<LO>(Ntotal*DofsPerNode+Nghost+1); Layerdofs = TLayerdofs.getRawPtr();
407 Teuchos::ArrayRCP<LO> TCol2Dof= Teuchos::arcp<LO>(Ntotal*DofsPerNode+Nghost+1); Col2Dof= TCol2Dof.getRawPtr();
409 RCP<Xpetra::Vector<LO,LO,GO,NO> > localdtemp = Xpetra::VectorFactory<LO,LO,GO,NO>::Build(Amat->getDomainMap());
410 RCP<Xpetra::Vector<LO,LO,GO,NO> > dtemp = Xpetra::VectorFactory<LO,LO,GO,NO>::Build(Amat->getColMap());
411 ArrayRCP<LO> valptr= localdtemp->getDataNonConst(0);
413 for (i = 0; i < Ntotal*DofsPerNode; i++)
414 valptr[i]= LayerId[i/DofsPerNode];
416 RCP< const Import> importer;
417 importer = Amat->getCrsGraph()->getImporter();
418 if (importer == Teuchos::null) {
419 importer = ImportFactory::Build(Amat->getDomainMap(), Amat->getColMap());
421 dtemp->doImport(*localdtemp, *(importer), Xpetra::INSERT);
423 valptr= dtemp->getDataNonConst(0);
424 for (i = 0; i < Ntotal*DofsPerNode+Nghost; i++) Layerdofs[i]= valptr[i];
425 valptr= localdtemp->getDataNonConst(0);
426 for (i = 0; i < Ntotal*DofsPerNode; i++) valptr[i]= i%DofsPerNode;
427 dtemp->doImport(*localdtemp, *(importer), Xpetra::INSERT);
428 valptr= dtemp->getDataNonConst(0);
429 for (i = 0; i < Ntotal*DofsPerNode+Nghost; i++) Col2Dof[i]= valptr[i];
432 NLayers = LayerId[0];
433 NVertLines= VertLineId[0];
435 else { NLayers = -1; NVertLines = -1; }
437 for (i = 1; i < Ntotal; i++) {
438 if ( VertLineId[i] > NVertLines ) NVertLines = VertLineId[i];
439 if ( LayerId[i] > NLayers ) NLayers = LayerId[i];
449 Teuchos::ArrayRCP<LO> TInvLineLayer= Teuchos::arcp<LO>(1+NVertLines*NLayers); InvLineLayer = TInvLineLayer.getRawPtr();
450 for (i=0; i < Ntotal; i++) {
451 InvLineLayer[ VertLineId[i]+1+LayerId[i]*NVertLines ] = i;
458 Teuchos::ArrayRCP<LO> TCptLayers = Teuchos::arcp<LO>(nz+1); CptLayers = TCptLayers.getRawPtr();
459 NCLayers = FindCpts(nz,CoarsenRate,0, &CptLayers);
468 if (NCLayers < 2) MaxStencilSize = nz;
469 else MaxStencilSize = CptLayers[2];
471 for (i = 3; i <= NCLayers; i++) {
472 if (MaxStencilSize < CptLayers[i]- CptLayers[i-2])
473 MaxStencilSize = CptLayers[i]- CptLayers[i-2];
476 if (MaxStencilSize < nz - CptLayers[NCLayers-1]+1)
477 MaxStencilSize = nz - CptLayers[NCLayers-1]+1;
487 Teuchos::ArrayRCP<LO> TSub2FullMap= Teuchos::arcp<LO>((MaxStencilSize+1)*DofsPerNode); Sub2FullMap= TSub2FullMap.getRawPtr();
488 Teuchos::ArrayRCP<SC> TBandSol= Teuchos::arcp<SC>((MaxStencilSize+1)*DofsPerNode*DofsPerNode); BandSol = TBandSol.getRawPtr();
492 KL = 2*DofsPerNode-1;
493 KU = 2*DofsPerNode-1;
497 Teuchos::ArrayRCP<SC> TBandMat= Teuchos::arcp<SC>(LDAB*MaxStencilSize*DofsPerNode+1); BandMat = TBandMat.getRawPtr();
498 Teuchos::ArrayRCP<LO> TIPIV= Teuchos::arcp<LO>((MaxStencilSize+1)*DofsPerNode); IPIV = TIPIV.getRawPtr();
507 Ndofs = DofsPerNode*Ntotal;
508 MaxNnz = 2*DofsPerNode*Ndofs;
510 RCP<const Map> rowMap = Amat->getRowMap();
511 Xpetra::global_size_t GNdofs= rowMap->getGlobalNumElements();
513 std::vector<size_t> stridingInfo_;
514 stridingInfo_.push_back(DofsPerNode);
516 Xpetra::global_size_t itemp = GNdofs/nz;
517 coarseMap = StridedMapFactory::Build(rowMap->lib(),
519 NCLayers*NVertLines*DofsPerNode,
530 P = rcp(
new CrsMatrixWrap(rowMap, coarseMap , 0));
531 coarseNullspace = MultiVectorFactory::Build(coarseMap, fineNullspace->getNumVectors());
534 Teuchos::ArrayRCP<SC> TPvals= Teuchos::arcp<SC>(1+MaxNnz); Pvals= TPvals.getRawPtr();
535 Teuchos::ArrayRCP<size_t> TPptr = Teuchos::arcp<size_t>(DofsPerNode*(2+Ntotal)); Pptr = TPptr.getRawPtr();
536 Teuchos::ArrayRCP<LO> TPcols= Teuchos::arcp<LO>(1+MaxNnz); Pcols= TPcols.getRawPtr();
538 Pptr[0] = 0; Pptr[1] = 0;
548 for (i = 1; i <= MaxNnz; i++) Pcols[i] = -1;
550 for (i = 1; i <= DofsPerNode*Ntotal+1; i++) {
552 count += (2*DofsPerNode);
564 for (MyLine=1; MyLine <= NVertLines; MyLine += 1) {
565 for (iii=1; iii <= NCLayers; iii+= 1) {
567 MyLayer = CptLayers[iii];
580 if (iii != 1 ) StartLayer = CptLayers[iii-1]+1;
583 if (iii != NCLayers) NStencilNodes= CptLayers[iii+1]-StartLayer;
584 else NStencilNodes= NLayers - StartLayer + 1;
587 N = NStencilNodes*DofsPerNode;
594 for (i = 0; i < NStencilNodes*DofsPerNode*DofsPerNode; i++)
596 for (i = 0; i < LDAB*N; i++) BandMat[ i] = 0.0;
603 for (node_k=1; node_k <= NStencilNodes ; node_k++) {
609 BlkRow = InvLineLayer[MyLine+(StartLayer+node_k-2)*NVertLines]+1;
610 Sub2FullMap[node_k] = BlkRow;
623 if (StartLayer+node_k-1 != MyLayer) {
624 for (
int dof_i=0; dof_i < DofsPerNode; dof_i++) {
626 j = (BlkRow-1)*DofsPerNode+dof_i;
627 ArrayView<const LO> AAcols;
628 ArrayView<const SC> AAvals;
629 Amat->getLocalRowView(j, AAcols, AAvals);
630 const int *Acols = AAcols.getRawPtr();
631 const SC *Avals = AAvals.getRawPtr();
632 RowLeng = AAvals.size();
634 TEUCHOS_TEST_FOR_EXCEPTION(RowLeng >= MaxNnzPerRow,
Exceptions::RuntimeError,
"MakeSemiCoarsenP: recompile with larger Max(HorNeighborNodes)\n");
636 for (i = 0; i < RowLeng; i++) {
637 LayDiff[i] = Layerdofs[Acols[i]]-StartLayer-node_k+2;
645 PtRow = (node_k-1)*DofsPerNode+dof_i+1;
646 for (dof_j=0; dof_j < DofsPerNode; dof_j++) {
647 PtCol = (node_k-1)*DofsPerNode+dof_j + 1;
651 for (i = 0; i < RowLeng; i++) {
652 if ((LayDiff[i] == 0) && (Col2Dof[Acols[i]] == dof_j))
655 index = LDAB*(PtCol-1)+KLU+PtRow-PtCol;
656 BandMat[index] = TheSum;
657 if (node_k != NStencilNodes) {
661 for (i = 0; i < RowLeng; i++) {
662 if ((LayDiff[i] == 1) &&(Col2Dof[Acols[i]]== dof_j))
665 j = PtCol+DofsPerNode;
666 index=LDAB*(j-1)+KLU+PtRow-j;
667 BandMat[index] = TheSum;
673 for (i = 0; i < RowLeng; i++) {
674 if ((LayDiff[i]== -1) &&(Col2Dof[Acols[i]]== dof_j))
677 j = PtCol-DofsPerNode;
678 index=LDAB*(j-1)+KLU+PtRow-j;
679 BandMat[index] = TheSum;
689 for (
int k = 0; k < static_cast<int>(fineNullspace->getNumVectors()); k++) {
690 Teuchos::ArrayRCP<SC> OneCNull = coarseNullspace->getDataNonConst(k);
691 Teuchos::ArrayRCP<SC> OneFNull = fineNullspace->getDataNonConst(k);
692 for (
int dof_i = 0; dof_i < DofsPerNode; dof_i++) {
693 OneCNull[( col-1)*DofsPerNode+dof_i] = OneFNull[ (BlkRow-1)*DofsPerNode+dof_i];
697 for (
int dof_i = 0; dof_i < DofsPerNode; dof_i++) {
700 PtRow = (node_k-1)*DofsPerNode+dof_i+1;
701 index = LDAB*(PtRow-1)+KLU;
702 BandMat[index] = 1.0;
703 BandSol[(dof_i)*DofsPerNode*NStencilNodes+PtRow-1] = 1.;
710 lapack.GBTRF( N, N, KL, KU, BandMat, LDAB, IPIV, &INFO);
714 lapack.GBTRS(notrans[0], N, KL, KU, NRHS, BandMat, LDAB, IPIV,
719 for (dof_j=0; dof_j < DofsPerNode; dof_j++) {
720 for (
int dof_i=0; dof_i < DofsPerNode; dof_i++) {
721 for (i =1; i <= NStencilNodes ; i++) {
722 index = (Sub2FullMap[i]-1)*DofsPerNode+dof_i+1;
724 Pcols[loc] = (col-1)*DofsPerNode+dof_j+1;
725 Pvals[loc] = BandSol[dof_j*DofsPerNode*NStencilNodes +
726 (i-1)*DofsPerNode + dof_i];
727 Pptr[index]= Pptr[index] + 1;
743 for (
size_t ii=1; ii <= Pptr[Ntotal*DofsPerNode]-1; ii++) {
744 if (ii == Pptr[CurRow]) {
745 Pptr[CurRow] = LastGuy;
747 while (ii > Pptr[CurRow]) {
748 Pptr[CurRow] = LastGuy;
752 if (Pcols[ii] != -1) {
753 Pcols[NewPtr-1] = Pcols[ii]-1;
754 Pvals[NewPtr-1] = Pvals[ii];
759 for (i = CurRow; i <= Ntotal*DofsPerNode; i++) Pptr[CurRow] = LastGuy;
764 for (i=-Ntotal*DofsPerNode+1; i>= 2 ; i--) {
765 Pptr[i-1] = Pptr[i-2];
769 ArrayRCP<size_t> rcpRowPtr;
770 ArrayRCP<LO> rcpColumns;
771 ArrayRCP<SC> rcpValues;
773 RCP<CrsMatrix> PCrs = rcp_dynamic_cast<CrsMatrixWrap>(P)->getCrsMatrix();
774 LO nnz = Pptr[Ndofs];
775 PCrs->allocateAllValues(nnz, rcpRowPtr, rcpColumns, rcpValues);
777 ArrayView<size_t> rowPtr = rcpRowPtr();
778 ArrayView<LO> columns = rcpColumns();
779 ArrayView<SC> values = rcpValues();
783 for (LO ii = 0; ii <= Ndofs; ii++) rowPtr[ii] = Pptr[ii];
784 for (LO ii = 0; ii < nnz; ii++) columns[ii] = Pcols[ii];
785 for (LO ii = 0; ii < nnz; ii++) values[ii] = Pvals[ii];
786 PCrs->setAllValues(rcpRowPtr, rcpColumns, rcpValues);
787 PCrs->expertStaticFillComplete(coarseMap, Amat->getDomainMap());
790 return NCLayers*NVertLines*DofsPerNode;
794 #define MUELU_SEMICOARSENPFACTORY_SHORT 795 #endif // MUELU_SEMICOARSENPFACTORY_DEF_HPP MueLu::DefaultLocalOrdinal LocalOrdinal
LO FindCpts(LO const PtsPerLine, LO const CoarsenRate, LO const Thin, LO **LayerCpts) const
Timer to be used in factories. Similar to Monitor but with additional timers.
void BuildP(Level &fineLevel, Level &coarseLevel) const
Abstract Build method.
Namespace for MueLu classes and methods.
static const NoFactory * get()
int GetLevelID() const
Return level number.
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.
Class that holds all level-specific information.
void Set(const std::string &ename, const T &entry, const FactoryBase *factory=NoFactory::get())
#define SET_VALID_ENTRY(name)
void Build(Level &fineLevel, Level &coarseLevel) const
Build an object with this factory.
#define MaxHorNeighborNodes
const RCP< const FactoryManagerBase > GetFactoryManager()
returns the current factory manager
bool IsAvailable(const std::string &ename, const FactoryBase *factory=NoFactory::get()) const
Test whether a need's value has been saved.
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()
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) const
void DeclareInput(Level &fineLevel, Level &coarseLevel) const
Input.