46 #ifndef MUELU_IPCFACTORY_DEF_HPP
47 #define MUELU_IPCFACTORY_DEF_HPP
60 #include "MueLu_PerfUtils.hpp"
62 #include "MueLu_Utilities.hpp"
72 #include "Intrepid2_HGRAD_LINE_C1_FEM.hpp"
73 #include "Intrepid2_HGRAD_LINE_Cn_FEM.hpp"
78 #include "Intrepid2_HGRAD_QUAD_C1_FEM.hpp"
80 #include "Intrepid2_HGRAD_QUAD_Cn_FEM.hpp"
95 #define MUELU_LEVEL_SET_IF_REQUESTED_OR_KEPT(level,ename,entry) \
96 {if (level.IsRequested(ename,this) || level.GetKeepFlag(ename,this) != 0) this->Set(level,ename,entry);}
104 namespace MueLuIntrepid {
105 inline std::string
tolower(
const std::string & str) {
106 std::string data(str);
107 std::transform(data.begin(), data.end(), data.begin(), [](
unsigned char c) {
return std::tolower(c); });
113 template<
class Basis,
class LOFieldContainer,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
115 std::vector<std::vector<LocalOrdinal> > &seeds,
129 shards::CellTopology cellTopo = basis->getBaseCellTopology();
130 int spaceDim = cellTopo.getDimension();
132 seeds.resize(spaceDim + 1);
133 typedef GlobalOrdinal
GO;
134 typedef LocalOrdinal
LO;
140 std::vector<std::set<LocalOrdinal>> seedSets(spaceDim+1);
142 int numCells = elementToNodeMap.extent(0);
143 for (
int cellOrdinal=0; cellOrdinal<numCells; cellOrdinal++)
145 for (
int d=0; d<=spaceDim; d++)
147 int subcellCount = cellTopo.getSubcellCount(d);
148 for (
int subcord=0; subcord<subcellCount; subcord++)
150 int dofCount = basis->getDofCount(d,subcord);
151 if (dofCount == 0)
continue;
153 GO leastGlobalDofOrdinal = go_invalid;
154 LO LID_leastGlobalDofOrdinal = lo_invalid;
155 for (
int basisOrdinalOrdinal=0; basisOrdinalOrdinal<dofCount; basisOrdinalOrdinal++)
157 int basisOrdinal = basis->getDofOrdinal(d,subcord,basisOrdinalOrdinal);
158 int colLID = elementToNodeMap(cellOrdinal,basisOrdinal);
163 if (rowLID != lo_invalid)
165 if ((leastGlobalDofOrdinal == go_invalid) || (colGID < leastGlobalDofOrdinal))
168 leastGlobalDofOrdinal = colGID;
169 LID_leastGlobalDofOrdinal = rowLID;
174 if (leastGlobalDofOrdinal != go_invalid)
176 seedSets[d].insert(LID_leastGlobalDofOrdinal);
181 for (
int d=0; d<=spaceDim;d++)
183 seeds[d] = std::vector<LocalOrdinal>(seedSets[d].begin(),seedSets[d].end());
194 template<
class Scalar,
class KokkosExecutionSpace>
199 string myerror(
"IntrepidBasisFactory: cannot parse string name '"+name+
"'");
204 size_t pos1 = name.find_first_of(
" _");
205 if(pos1==0)
throw std::runtime_error(myerror);
206 string deriv =
tolower(name.substr(0,pos1));
207 if(deriv!=
"hgrad" && deriv!=
"hcurl" && deriv!=
"hdiv")
throw std::runtime_error(myerror);
211 size_t pos2 = name.find_first_of(
" _",pos1);
212 if(pos2==0)
throw std::runtime_error(myerror);
213 string el =
tolower(name.substr(pos1,pos2-pos1));
214 if(el!=
"hex" && el!=
"line" && el!=
"poly" && el!=
"pyr" && el!=
"quad" && el!=
"tet" && el!=
"tri" && el!=
"wedge")
throw std::runtime_error(myerror);
218 string poly =
tolower(name.substr(pos2,1));
219 if(poly!=
"c" && poly!=
"i")
throw std::runtime_error(myerror);
223 degree=std::stoi(name.substr(pos2,1));
224 if(degree<=0)
throw std::runtime_error(myerror);
227 if(deriv==
"hgrad" && el==
"quad" && poly==
"c"){
228 if(degree==1)
return rcp(
new Intrepid2::Basis_HGRAD_QUAD_C1_FEM<KokkosExecutionSpace,Scalar,Scalar>());
229 else return rcp(
new Intrepid2::Basis_HGRAD_QUAD_Cn_FEM<KokkosExecutionSpace,Scalar,Scalar>(degree,Intrepid2::POINTTYPE_EQUISPACED));
231 else if(deriv==
"hgrad" && el==
"line" && poly==
"c"){
232 if(degree==1)
return rcp(
new Intrepid2::Basis_HGRAD_LINE_C1_FEM<KokkosExecutionSpace,Scalar,Scalar>());
233 else return rcp(
new Intrepid2::Basis_HGRAD_LINE_Cn_FEM<KokkosExecutionSpace,Scalar,Scalar>(degree,Intrepid2::POINTTYPE_EQUISPACED));
237 throw std::runtime_error(myerror);
248 template<
class Scalar,
class KokkosDeviceType>
250 std::vector<size_t> & lo_node_in_hi,
251 Kokkos::DynRankView<Scalar,KokkosDeviceType> & hi_DofCoords) {
253 typedef typename KokkosDeviceType::execution_space KokkosExecutionSpace;
255 size_t degree = hi_basis->getDegree();
256 lo_node_in_hi.resize(0);
258 if(!rcp_dynamic_cast<Intrepid2::Basis_HGRAD_QUAD_Cn_FEM<KokkosExecutionSpace,Scalar,Scalar> >(hi_basis).
is_null()) {
260 lo_node_in_hi.insert(lo_node_in_hi.end(),{0,degree, (degree+1)*(degree+1)-1, degree*(degree+1)});
262 else if(!rcp_dynamic_cast<Intrepid2::Basis_HGRAD_LINE_Cn_FEM<KokkosExecutionSpace,Scalar,Scalar> >(hi_basis).is_null()) {
264 lo_node_in_hi.insert(lo_node_in_hi.end(),{0,degree});
267 throw std::runtime_error(
"IntrepidPCoarsenFactory: Unknown element type");
270 Kokkos::resize(hi_DofCoords,hi_basis->getCardinality(),hi_basis->getBaseCellTopology().getDimension());
271 hi_basis->getDofCoords(hi_DofCoords);
283 template<
class LocalOrdinal,
class GlobalOrdinal,
class Node,
class LOFieldContainer>
286 LOFieldContainer & lo_elemToHiRepresentativeNode) {
287 typedef GlobalOrdinal
GO;
292 size_t numElem = hi_elemToNode.extent(0);
293 size_t lo_nperel = candidates.size();
297 for(
size_t i=0; i<numElem; i++)
298 for(
size_t j=0; j<lo_nperel; j++) {
299 if(candidates[j].size() == 1)
300 lo_elemToHiRepresentativeNode(i,j) = hi_elemToNode(i,candidates[j][0]);
303 std::vector<GO> GID(candidates[j].size());
304 for(
size_t k=0; k<(size_t)candidates[j].size(); k++)
305 GID[k] = hi_columnMap->getGlobalElement(hi_elemToNode(i,candidates[j][k]));
308 size_t which = std::distance(GID.begin(),std::min_element(GID.begin(),GID.end()));
311 lo_elemToHiRepresentativeNode(i,j) = hi_elemToNode(i,candidates[j][which]);
326 template <
class LocalOrdinal,
class LOFieldContainer>
328 const std::vector<bool> & hi_nodeIsOwned,
329 const LOFieldContainer & lo_elemToHiRepresentativeNode,
330 LOFieldContainer & lo_elemToNode,
331 std::vector<bool> & lo_nodeIsOwned,
332 std::vector<LocalOrdinal> & hi_to_lo_map,
333 int & lo_numOwnedNodes) {
334 typedef LocalOrdinal
LO;
337 size_t numElem = hi_elemToNode.extent(0);
338 size_t hi_numNodes = hi_nodeIsOwned.size();
339 size_t lo_nperel = lo_elemToHiRepresentativeNode.extent(1);
343 std::vector<bool> is_low_order(hi_numNodes,
false);
344 for(
size_t i=0; i<numElem; i++)
345 for(
size_t j=0; j<lo_nperel; j++) {
346 LO
id = lo_elemToHiRepresentativeNode(i,j);
347 is_low_order[id] =
true;
352 size_t lo_numNodes=0;
355 for(
size_t i=0; i<hi_numNodes; i++)
356 if(is_low_order[i]) {
357 hi_to_lo_map[i] = lo_numNodes;
359 if(hi_nodeIsOwned[i]) lo_numOwnedNodes++;
363 lo_nodeIsOwned.resize(lo_numNodes,
false);
364 for(
size_t i=0; i<hi_numNodes; i++) {
365 if(is_low_order[i] && hi_nodeIsOwned[i])
366 lo_nodeIsOwned[hi_to_lo_map[i]]=
true;
370 for(
size_t i=0; i<numElem; i++)
371 for(
size_t j=0; j<lo_nperel; j++)
372 lo_elemToNode(i,j) = hi_to_lo_map[lo_elemToHiRepresentativeNode(i,j)];
377 bool map_ordering_test_passed=
true;
378 for(
size_t i=0; i<lo_numNodes-1; i++)
379 if(!lo_nodeIsOwned[i] && lo_nodeIsOwned[i+1])
380 map_ordering_test_passed=
false;
382 if(!map_ordering_test_passed)
383 throw std::runtime_error(
"MueLu::MueLuIntrepid::BuildLoElemToNodeViaRepresentatives failed map ordering test");
399 template <
class LocalOrdinal,
class LOFieldContainer>
401 const std::vector<bool> & hi_nodeIsOwned,
402 const std::vector<size_t> & lo_node_in_hi,
404 LOFieldContainer & lo_elemToNode,
405 std::vector<bool> & lo_nodeIsOwned,
406 std::vector<LocalOrdinal> & hi_to_lo_map,
407 int & lo_numOwnedNodes) {
408 typedef LocalOrdinal
LO;
413 size_t numElem = hi_elemToNode.extent(0);
414 size_t hi_numNodes = hi_nodeIsOwned.size();
416 size_t lo_nperel = lo_node_in_hi.size();
420 std::vector<bool> is_low_order(hi_numNodes,
false);
421 for(
size_t i=0; i<numElem; i++)
422 for(
size_t j=0; j<lo_nperel; j++) {
423 LO lid = hi_elemToNode(i,lo_node_in_hi[j]);
426 if(hi_isDirichlet[lid])
427 lo_elemToNode(i,j) = LOINVALID;
429 lo_elemToNode(i,j) = lid;
430 is_low_order[hi_elemToNode(i,lo_node_in_hi[j])] =
true;
436 size_t lo_numNodes=0;
439 for(
size_t i=0; i<hi_numNodes; i++)
440 if(is_low_order[i]) {
441 hi_to_lo_map[i] = lo_numNodes;
443 if(hi_nodeIsOwned[i]) lo_numOwnedNodes++;
447 lo_nodeIsOwned.resize(lo_numNodes,
false);
448 for(
size_t i=0; i<hi_numNodes; i++) {
449 if(is_low_order[i] && hi_nodeIsOwned[i])
450 lo_nodeIsOwned[hi_to_lo_map[i]]=
true;
454 for(
size_t i=0; i<numElem; i++)
455 for(
size_t j=0; j<lo_nperel; j++) {
456 if(lo_elemToNode(i,j) != LOINVALID)
457 lo_elemToNode(i,j) = hi_to_lo_map[lo_elemToNode(i,j)];
462 bool map_ordering_test_passed=
true;
463 for(
size_t i=0; i<lo_numNodes-1; i++)
464 if(!lo_nodeIsOwned[i] && lo_nodeIsOwned[i+1])
465 map_ordering_test_passed=
false;
467 if(!map_ordering_test_passed)
468 throw std::runtime_error(
"MueLu::MueLuIntrepid::BuildLoElemToNode failed map ordering test");
481 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
482 void GenerateColMapFromImport(
const Xpetra::Import<LocalOrdinal,GlobalOrdinal, Node> & hi_importer,
const std::vector<LocalOrdinal> &hi_to_lo_map,
const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node> & lo_domainMap,
const size_t & lo_columnMapLength,
RCP<
const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > & lo_columnMap) {
483 typedef LocalOrdinal
LO;
484 typedef GlobalOrdinal
GO;
501 for(
size_t i=0; i<hi_domainMap->getNodeNumElements(); i++) {
502 if(hi_to_lo_map[i]!=lo_invalid) dvec_data[i] = lo_domainMap.
getGlobalElement(hi_to_lo_map[i]);
503 else dvec_data[i] = go_invalid;
512 Array<GO> lo_col_data(lo_columnMapLength);
514 for(
size_t i=0,idx=0; i<hi_columnMap->getNodeNumElements(); i++) {
515 if(hi_to_lo_map[i]!=lo_invalid) {
516 lo_col_data[idx] = cvec_data[i];
532 template<
class Basis,
class SCFieldContainer>
533 void GenerateRepresentativeBasisNodes(
const Basis & basis,
const SCFieldContainer & ReferenceNodeLocations,
const double threshold,std::vector<std::vector<size_t> > & representative_node_candidates) {
534 typedef SCFieldContainer FC;
535 typedef typename FC::data_type
SC;
538 size_t numFieldsHi = ReferenceNodeLocations.extent(0);
540 size_t numFieldsLo = basis.getCardinality();
542 FC LoValues(
"LoValues",numFieldsLo,numFieldsHi);
544 basis.getValues(LoValues, ReferenceNodeLocations , Intrepid2::OPERATOR_VALUE);
548 printf(
"** LoValues[%d,%d] **\n",(
int)numFieldsLo,(
int)numFieldsHi);
549 for(
size_t i=0; i<numFieldsLo; i++) {
550 for(
size_t j=0; j<numFieldsHi; j++)
551 printf(
"%6.4e ",LoValues(i,j));
554 printf(
"**************\n");fflush(stdout);
557 representative_node_candidates.resize(numFieldsLo);
558 for(
size_t i=0; i<numFieldsLo; i++) {
561 for(
size_t j=0; j<numFieldsHi; j++)
565 for(
size_t j=0; j<numFieldsHi; j++) {
567 representative_node_candidates[i].push_back(j);
572 for(
size_t i=0; i<numFieldsLo; i++)
573 if(!representative_node_candidates[i].size())
574 throw std::runtime_error(
"ERROR: GenerateRepresentativeBasisNodes: No candidates found!");
587 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
589 const std::vector<bool> & hi_nodeIsOwned,
591 const std::vector<size_t> &lo_node_in_hi,
592 const Basis &lo_basis,
593 const std::vector<LocalOrdinal> & hi_to_lo_map,
600 size_t numFieldsHi = hi_elemToNode.extent(1);
601 size_t numFieldsLo = lo_basis.getCardinality();
603 FC LoValues_at_HiDofs(
"LoValues_at_HiDofs",numFieldsLo,numFieldsHi);
604 lo_basis.getValues(LoValues_at_HiDofs, hi_DofCoords, Intrepid2::OPERATOR_VALUE);
611 P =
rcp(
new CrsMatrixWrap(hi_map,lo_colMap,0));
612 RCP<CrsMatrix> Pcrs = rcp_dynamic_cast<CrsMatrixWrap>(P)->getCrsMatrix();
615 size_t Nelem=hi_elemToNode.extent(0);
616 std::vector<bool> touched(hi_map->getNodeNumElements(),
false);
619 for(
size_t i=0; i<Nelem; i++) {
620 for(
size_t j=0; j<numFieldsHi; j++) {
621 LO row_lid = hi_elemToNode(i,j);
622 GO row_gid = hi_map->getGlobalElement(row_lid);
623 if(hi_nodeIsOwned[row_lid] && !touched[row_lid]) {
624 for(
size_t k=0; k<numFieldsLo; k++) {
626 LO col_lid = hi_to_lo_map[hi_elemToNode(i,lo_node_in_hi[k])];
627 if(col_lid==LOINVALID)
continue;
629 col_gid[0] = {lo_colMap->getGlobalElement(col_lid)};
630 val[0] = LoValues_at_HiDofs(k,j);
634 P->insertGlobalValues(row_gid,col_gid(),val());
636 touched[row_lid]=
true;
640 P->fillComplete(lo_domainMap,hi_map);
644 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
646 const std::vector<bool> & hi_nodeIsOwned,
649 const Basis &lo_basis,
650 const std::vector<LocalOrdinal> & hi_to_lo_map,
657 size_t numFieldsHi = hi_elemToNode.extent(1);
658 size_t numFieldsLo = lo_basis.getCardinality();
659 FC LoValues_at_HiDofs(
"LoValues_at_HiDofs",numFieldsLo,numFieldsHi);
660 lo_basis.getValues(LoValues_at_HiDofs, hi_DofCoords, Intrepid2::OPERATOR_VALUE);
667 P =
rcp(
new CrsMatrixWrap(hi_map,lo_colMap,0));
668 RCP<CrsMatrix> Pcrs = rcp_dynamic_cast<CrsMatrixWrap>(P)->getCrsMatrix();
671 size_t Nelem=hi_elemToNode.extent(0);
672 std::vector<bool> touched(hi_map->getNodeNumElements(),
false);
675 for(
size_t i=0; i<Nelem; i++) {
676 for(
size_t j=0; j<numFieldsHi; j++) {
677 LO row_lid = hi_elemToNode(i,j);
678 GO row_gid = hi_map->getGlobalElement(row_lid);
679 if(hi_nodeIsOwned[row_lid] && !touched[row_lid]) {
680 for(
size_t k=0; k<numFieldsLo; k++) {
682 LO col_lid = hi_to_lo_map[lo_elemToHiRepresentativeNode(i,k)];
683 col_gid[0] = {lo_colMap->getGlobalElement(col_lid)};
684 val[0] = LoValues_at_HiDofs(k,j);
688 P->insertGlobalValues(row_gid,col_gid(),val());
690 touched[row_lid]=
true;
694 P->fillComplete(lo_domainMap,hi_map);
698 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
702 #define SET_VALID_ENTRY(name) validParamList->setEntry(name, MasterList::getEntry(name))
705 #undef SET_VALID_ENTRY
707 validParamList->set<
RCP<const FactoryBase> >(
"A", Teuchos::null,
"Generating factory of the matrix A used during the prolongator smoothing process");
709 validParamList->set<
RCP<const FactoryBase> >(
"Nullspace", Teuchos::null,
"Generating factory of the nullspace");
710 validParamList->set<
RCP<const FactoryBase> >(
"pcoarsen: element to node map", Teuchos::null,
"Generating factory of the element to node map");
711 return validParamList;
715 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
717 Input(fineLevel,
"A");
718 Input(fineLevel,
"pcoarsen: element to node map");
719 Input(fineLevel,
"Nullspace");
723 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
725 return BuildP(fineLevel, coarseLevel);
729 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
733 const std::string prefix =
"MueLu::IntrepidPCoarsenFactory(" + levelIDs +
"): ";
736 typedef Kokkos::DynRankView<LocalOrdinal,typename Node::device_type> FCi;
737 typedef Kokkos::DynRankView<double,typename Node::device_type> FC;
740 RCP<Matrix> A = Get< RCP<Matrix> >(fineLevel,
"A");
741 RCP<MultiVector> fineNullspace = Get< RCP<MultiVector> >(fineLevel,
"Nullspace");
745 if (restrictionMode_) {
751 std::vector<LocalOrdinal> A_dirichletRows;
759 if (coarseLevel.
IsAvailable(
"AP reuse data",
this)) {
760 GetOStream(static_cast<MsgType>(
Runtime0 |
Test)) <<
"Reusing previous AP data" << std::endl;
764 if (APparams->isParameter(
"graph"))
774 int lo_degree, hi_degree;
775 RCP<Basis> hi_basis = MueLuIntrepid::BasisFactory<double,typename Node::device_type::execution_space>(pL.
get<std::string>(
"pcoarsen: hi basis"),hi_degree);
776 RCP<Basis> lo_basis = MueLuIntrepid::BasisFactory<double,typename Node::device_type::execution_space>(pL.
get<std::string>(
"pcoarsen: lo basis"),lo_degree);
779 GetOStream(
Statistics1) <<
"P-Coarsening from basis "<<pL.
get<std::string>(
"pcoarsen: hi basis")<<
" to "<<pL.
get<std::string>(
"pcoarsen: lo basis") <<std::endl;
783 const Teuchos::RCP<FCi> Pn_elemToNode = Get<Teuchos::RCP<FCi> >(fineLevel,
"pcoarsen: element to node map");
793 int NumProc = rowMap->getComm()->getSize();
794 assert(rowMap->isSameAs(*domainMap));
795 std::vector<bool> Pn_nodeIsOwned(colMap->getNodeNumElements(),
false);
796 LO num_owned_rows = 0;
797 for(
size_t i=0; i<rowMap->getNodeNumElements(); i++) {
798 if(rowMap->getGlobalElement(i) == colMap->getGlobalElement(i)) {
799 Pn_nodeIsOwned[i] =
true;
808 std::vector<bool> P1_nodeIsOwned;
809 int P1_numOwnedNodes;
810 std::vector<LO> hi_to_lo_map;
813 std::vector<size_t> lo_node_in_hi;
816 FCi lo_elemToHiRepresentativeNode;
823 printf(
"[%d] isDirichletRow = ",A->getRowMap()->getComm()->getRank());
824 for(
size_t i=0;i<hi_isDirichletRow->getMap()->getNodeNumElements(); i++)
825 printf(
"%d ",hi_isDirichletRow->getData(0)[i]);
827 printf(
"[%d] isDirichletCol = ",A->getRowMap()->getComm()->getRank());
828 for(
size_t i=0;i<hi_isDirichletCol->getMap()->getNodeNumElements(); i++)
829 printf(
"%d ",hi_isDirichletCol->getData(0)[i]);
840 MueLuIntrepid::BuildLoElemToNode(*Pn_elemToNode,Pn_nodeIsOwned,lo_node_in_hi,hi_isDirichletCol->getData(0),*P1_elemToNode,P1_nodeIsOwned,hi_to_lo_map,P1_numOwnedNodes);
841 assert(hi_to_lo_map.size() == colMap->getNodeNumElements());
845 double threshold = 1e-10;
846 std::vector<std::vector<size_t> > candidates;
847 Kokkos::resize(hi_DofCoords,hi_basis->getCardinality(),hi_basis->getBaseCellTopology().getDimension());
848 hi_basis->getDofCoords(hi_DofCoords);
850 MueLu::MueLuIntrepid::GenerateRepresentativeBasisNodes<Basis,FC>(*lo_basis,hi_DofCoords,threshold,candidates);
866 if(NumProc==1) P1_colMap = P1_domainMap;
867 else MueLuIntrepid::GenerateColMapFromImport<LO,GO,NO>(*Acrs.getCrsGraph()->getImporter(),hi_to_lo_map,*P1_domainMap,P1_nodeIsOwned.size(),P1_colMap);
872 GenerateLinearCoarsening_pn_kirby_to_p1(*Pn_elemToNode,Pn_nodeIsOwned,hi_DofCoords,lo_node_in_hi,*lo_basis,hi_to_lo_map,P1_colMap,P1_domainMap,A->getRowMap(),finalP);
874 GenerateLinearCoarsening_pn_kirby_to_pm(*Pn_elemToNode,Pn_nodeIsOwned,hi_DofCoords,lo_elemToHiRepresentativeNode,*lo_basis,hi_to_lo_map,P1_colMap,P1_domainMap,A->getRowMap(),finalP);
882 RCP<MultiVector> coarseNullspace = MultiVectorFactory::Build(P1_domainMap, fineNullspace->getNumVectors());
884 Set(coarseLevel,
"Nullspace", coarseNullspace);
887 if (!restrictionMode_) {
889 Set(coarseLevel,
"P", finalP);
891 APparams->set(
"graph", finalP);
896 params->set(
"printLoadBalancingInfo",
true);
897 params->set(
"printCommInfo",
true);
908 Set(coarseLevel,
"R", R);
912 params->set(
"printLoadBalancingInfo",
true);
913 params->set(
"printCommInfo",
true);
922 #endif // MUELU_IPCFACTORY_DEF_HPP
void GenerateLoNodeInHiViaGIDs(const std::vector< std::vector< size_t > > &candidates, const LOFieldContainer &hi_elemToNode, RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > &hi_columnMap, LOFieldContainer &lo_elemToHiRepresentativeNode)
virtual Teuchos::RCP< const Teuchos::Comm< int > > getComm() const =0
T & Get(const std::string &ename, const FactoryBase *factory=NoFactory::get())
Get data without decrementing associated storage counter (i.e., read-only access). Usage: Level->Get< RCP<Matrix> >("A", factory) if factory == NULL => use default factory.
Teuchos::RCP< Intrepid2::Basis< KokkosExecutionSpace, Scalar, Scalar > > BasisFactory(const std::string &name, int °ree)
bool is_null(const boost::shared_ptr< T > &p)
void FindGeometricSeedOrdinals(Teuchos::RCP< Basis > basis, const LOFieldContainer &elementToNodeMap, std::vector< std::vector< LocalOrdinal > > &seeds, const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &rowMap, const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &columnMap)
std::string toString(const T &what)
Little helper function to convert non-string types to strings.
static magnitudeType eps()
T & get(const std::string &name, T def_value)
void BuildP(Level &fineLevel, Level &coarseLevel) const
Abstract Build method.
Timer to be used in factories. Similar to Monitor but with additional timers.
One-liner description of what is happening.
virtual LocalOrdinal getLocalElement(GlobalOrdinal globalIndex) const =0
std::string tolower(const std::string &str)
virtual GlobalOrdinal getIndexBase() const =0
std::enable_if< std::is_same< typename Kokkos::View< T, P...>::array_layout, Kokkos::LayoutLeft >::value||std::is_same< typename Kokkos::View< T, P...>::array_layout, Kokkos::LayoutRight >::value >::type resize(Kokkos::View< T, P...> &v, const size_t n0=KOKKOS_IMPL_CTOR_DEFAULT_ARG, const size_t n1=KOKKOS_IMPL_CTOR_DEFAULT_ARG, const size_t n2=KOKKOS_IMPL_CTOR_DEFAULT_ARG, const size_t n3=KOKKOS_IMPL_CTOR_DEFAULT_ARG, const size_t n4=KOKKOS_IMPL_CTOR_DEFAULT_ARG, const size_t n5=KOKKOS_IMPL_CTOR_DEFAULT_ARG, const size_t n6=KOKKOS_IMPL_CTOR_DEFAULT_ARG, const size_t n7=KOKKOS_IMPL_CTOR_DEFAULT_ARG)
RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > getDomainMap() const
Print even more statistics.
virtual Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getTargetMap() const =0
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Class that holds all level-specific information.
Kokkos::DynRankView< double, typename Node::device_type > SCFieldContainer
static void FindDirichletRowsAndPropagateToCols(Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, Teuchos::RCP< Xpetra::Vector< int, LocalOrdinal, GlobalOrdinal, Node > > &isDirichletRow, Teuchos::RCP< Xpetra::Vector< int, LocalOrdinal, GlobalOrdinal, Node > > &isDirichletCol)
Timer to be used in factories. Similar to SubMonitor but adds a timer level by level.
void BuildLoElemToNode(const LOFieldContainer &hi_elemToNode, const std::vector< bool > &hi_nodeIsOwned, const std::vector< size_t > &lo_node_in_hi, const Teuchos::ArrayRCP< const int > &hi_isDirichlet, LOFieldContainer &lo_elemToNode, std::vector< bool > &lo_nodeIsOwned, std::vector< LocalOrdinal > &hi_to_lo_map, int &lo_numOwnedNodes)
void BuildLoElemToNodeViaRepresentatives(const LOFieldContainer &hi_elemToNode, const std::vector< bool > &hi_nodeIsOwned, const LOFieldContainer &lo_elemToHiRepresentativeNode, LOFieldContainer &lo_elemToNode, std::vector< bool > &lo_nodeIsOwned, std::vector< LocalOrdinal > &hi_to_lo_map, int &lo_numOwnedNodes)
#define MUELU_LEVEL_SET_IF_REQUESTED_OR_KEPT(level, ename, entry)
static void FindDirichletRows(Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, std::vector< LocalOrdinal > &dirichletRows, bool count_twos_as_dirichlet=false)
static RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Transpose(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, bool optimizeTranspose=false, const std::string &label=std::string(), const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.
void GenerateLinearCoarsening_pn_kirby_to_p1(const LOFieldContainer &hi_elemToNode, const std::vector< bool > &hi_nodeIsOwned, const SCFieldContainer &hi_DofCoords, const std::vector< size_t > &lo_node_in_hi, const Basis &lo_Basis, const std::vector< LocalOrdinal > &hi_to_lo_map, const Teuchos::RCP< const Map > &lo_colMap, const Teuchos::RCP< const Map > &lo_domainMap, const Teuchos::RCP< const Map > &hi_map, Teuchos::RCP< Matrix > &P) const
virtual GlobalOrdinal getGlobalElement(LocalOrdinal localIndex) const =0
void DeclareInput(Level &fineLevel, Level &coarseLevel) const
Input.
static std::string PrintMatrixInfo(const Matrix &A, const std::string &msgTag, RCP< const Teuchos::ParameterList > params=Teuchos::null)
static void ZeroDirichletRows(Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const std::vector< LocalOrdinal > &dirichletRows, Scalar replaceWith=Teuchos::ScalarTraits< Scalar >::zero())
static RCP< Vector > Build(const Teuchos::RCP< const Map > &map, bool zeroOut=true)
Kokkos::DynRankView< LocalOrdinal, typename Node::device_type > LOFieldContainer
Intrepid2::Basis< typename Node::device_type::execution_space, double, double > Basis
void GenerateColMapFromImport(const Xpetra::Import< LocalOrdinal, GlobalOrdinal, Node > &hi_importer, const std::vector< LocalOrdinal > &hi_to_lo_map, const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &lo_domainMap, const size_t &lo_columnMapLength, RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > &lo_columnMap)
#define TEUCHOS_UNREACHABLE_RETURN(dummyReturnVal)
void Build(Level &fineLevel, Level &coarseLevel) const
Build method.
void IntrepidGetP1NodeInHi(const Teuchos::RCP< Intrepid2::Basis< typename KokkosDeviceType::execution_space, Scalar, Scalar > > &hi_basis, std::vector< size_t > &lo_node_in_hi, Kokkos::DynRankView< Scalar, KokkosDeviceType > &hi_DofCoords)
void GenerateRepresentativeBasisNodes(const Basis &basis, const SCFieldContainer &ReferenceNodeLocations, const double threshold, std::vector< std::vector< size_t > > &representative_node_candidates)
int GetLevelID() const
Return level number.
void GenerateLinearCoarsening_pn_kirby_to_pm(const LOFieldContainer &hi_elemToNode, const std::vector< bool > &hi_nodeIsOwned, const SCFieldContainer &hi_DofCoords, const LOFieldContainer &lo_elemToHiRepresentativeNode, const Basis &lo_basis, const std::vector< LocalOrdinal > &hi_to_lo_map, const Teuchos::RCP< const Map > &lo_colMap, const Teuchos::RCP< const Map > &lo_domainMap, const Teuchos::RCP< const Map > &hi_map, Teuchos::RCP< Matrix > &P) const
static Teuchos::RCP< Map< LocalOrdinal, GlobalOrdinal, Node > > Build(UnderlyingLib lib, global_size_t numGlobalElements, GlobalOrdinal indexBase, const Teuchos::RCP< const Teuchos::Comm< int > > &comm, LocalGlobal lg=Xpetra::GloballyDistributed)
virtual Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getSourceMap() const =0
virtual UnderlyingLib lib() const =0
#define SET_VALID_ENTRY(name)
bool IsAvailable(const std::string &ename, const FactoryBase *factory=NoFactory::get()) const
Test whether a need's value has been saved.