43 #ifndef PANZER_DOF_MANAGER2_IMPL_HPP
44 #define PANZER_DOF_MANAGER2_IMPL_HPP
51 #include "PanzerDofMgr_config.hpp"
64 #include "Teuchos_ArrayView.hpp"
66 #include "Tpetra_Map.hpp"
67 #include "Tpetra_Export.hpp"
68 #include "Tpetra_Vector.hpp"
69 #include "Tpetra_MultiVector.hpp"
71 #include <unordered_set>
73 #ifdef PHX_KOKKOS_DEVICE_TYPE_CUDA
74 #define PANZER_DOFMGR_REQUIRE_CUDA
88 template <
typename LocalOrdinal,
typename GlobalOrdinal>
89 class HashTieBreak :
public Tpetra::Details::TieBreak<LocalOrdinal,GlobalOrdinal> {
93 HashTieBreak(
const unsigned int seed = (2654435761U))
96 virtual std::size_t selectedIndex(GlobalOrdinal GID,
97 const std::vector<std::pair<int,LocalOrdinal> > & pid_and_lid)
const
100 int intkey = (int) ((GID & 0x000000007fffffffLL) +
101 ((GID & 0x7fffffff80000000LL) >> 31));
102 return std::size_t((
seed_ ^ intkey) % pid_and_lid.size());
110 template <
typename LocalOrdinal,
typename GlobalOrdinal>
111 class GreedyTieBreak :
public Tpetra::Details::TieBreak<LocalOrdinal,GlobalOrdinal> {
116 virtual bool mayHaveSideEffects()
const {
120 virtual std::size_t selectedIndex(GlobalOrdinal ,
121 const std::vector<std::pair<int,LocalOrdinal> > & pid_and_lid)
const
124 const std::size_t numLids = pid_and_lid.size();
126 int minpid = pid_and_lid[0].first;
127 std::size_t minidx = 0;
128 for (idx = 0; idx < numLids; ++idx) {
129 if (pid_and_lid[idx].first < minpid) {
130 minpid = pid_and_lid[idx].first;
147 template <
typename LO,
typename GO>
149 : numFields_(0),buildConnectivityRun_(false),requireOrientations_(false), useTieBreak_(false), useNeighbors_(false)
152 template <
typename LO,
typename GO>
154 : numFields_(0),buildConnectivityRun_(false),requireOrientations_(false), useTieBreak_(false), useNeighbors_(false)
159 template <
typename LO,
typename GO>
163 "DOFManager::setConnManager: setConnManager cannot be called after "
164 "buildGlobalUnknowns has been called");
165 connMngr_ = connMngr;
167 connMngr_->getElementBlockIds(blockOrder_);
168 for (
size_t i = 0; i < blockOrder_.size(); ++i) {
169 blockNameToID_.insert(std::map<std::string,int>::value_type(blockOrder_[i],i));
172 blockToAssociatedFP_.resize(blockOrder_.size());
179 template <
typename LO,
typename GO>
184 "DOFManager::addField: addField cannot be called after "
185 "buildGlobalUnknowns has been called");
187 fieldPatterns_.push_back(pattern);
188 fieldTypes_.push_back(type);
189 fieldNameToAID_.insert(std::map<std::string,int>::value_type(str, numFields_));
192 fieldStringOrder_.push_back(str);
193 fieldAIDOrder_.push_back(numFields_);
195 for(std::size_t i=0;i<blockOrder_.size();i++) {
196 blockToAssociatedFP_[i].push_back(numFields_);
203 template <
typename LO,
typename GO>
208 "DOFManager::addField: addField cannot be called after "
209 "buildGlobalUnknowns has been called");
211 "DOFManager::addField: you must add a ConnManager before"
212 "you can associate a FP with a given block.")
215 while(!found && blocknum<blockOrder_.size()){
216 if(blockOrder_[blocknum]==blockID){
227 std::map<std::string,int>::const_iterator fpIter = fieldNameToAID_.find(str);
228 if(fpIter!=fieldNameToAID_.end())
232 fieldPatterns_.push_back(pattern);
233 fieldTypes_.push_back(type);
234 fieldNameToAID_.insert(std::map<std::string,int>::value_type(str, numFields_));
236 fieldStringOrder_.push_back(str);
237 fieldAIDOrder_.push_back(numFields_);
240 blockToAssociatedFP_[blocknum].push_back(numFields_);
245 blockToAssociatedFP_[blocknum].push_back(fpIter->second);
253 template <
typename LO,
typename GO>
256 std::map<std::string,int>::const_iterator fitr = fieldNameToAID_.find(name);
257 if(fitr==fieldNameToAID_.end())
258 return Teuchos::null;
260 if(fitr->second<
int(fieldPatterns_.size()))
261 return fieldPatterns_[fitr->second];
263 return Teuchos::null;
267 template <
typename LO,
typename GO>
270 std::map<std::string,int>::const_iterator fitr = fieldNameToAID_.find(fieldName);
271 if(fitr==fieldNameToAID_.end())
272 return Teuchos::null;
276 while(!found && blocknum<blockOrder_.size()){
277 if(blockOrder_[blocknum]==blockId){
285 return Teuchos::null;
287 std::vector<int>::const_iterator itr = std::find(blockToAssociatedFP_[blocknum].begin(),
288 blockToAssociatedFP_[blocknum].
end(),fitr->second);
289 if(itr!=blockToAssociatedFP_[blocknum].
end()) {
290 if(fitr->second<
int(fieldPatterns_.size()))
291 return fieldPatterns_[fitr->second];
294 return Teuchos::null;
302 template<
typename LO,
typename GO>
306 std::vector<GO>& indices)
const
316 template<
typename LO,
typename GO>
320 std::vector<GO>& indices)
const
330 template<
typename LO,
typename GO>
334 std::vector<GO>& indices)
const
337 indices.resize(owned_.size() + ghosted_.size());
338 for (
size_t i(0); i < owned_.size(); ++i)
339 indices[i] = owned_[i];
340 for (
size_t i(0); i < ghosted_.size(); ++i)
341 indices[owned_.size() + i] = ghosted_[i];
349 template<
typename LO,
typename GO>
354 return owned_.size();
362 template<
typename LO,
typename GO>
367 return ghosted_.size();
375 template<
typename LO,
typename GO>
380 return owned_.size() + ghosted_.size();
384 template <
typename LO,
typename GO>
390 template <
typename LO,
typename GO>
393 TEUCHOS_TEST_FOR_EXCEPTION(!buildConnectivityRun_,std::logic_error,
"DOFManager::getGIDFieldOffsets: cannot be called before "
394 "buildGlobalUnknowns has been called");
395 std::map<std::string,int>::const_iterator bitr = blockNameToID_.find(blockID);
396 if(bitr==blockNameToID_.end())
398 int bid=bitr->second;
399 if(fa_fps_[bid]!=Teuchos::null)
400 return fa_fps_[bid]->localOffsets(fieldNum);
402 static const std::vector<int> empty;
406 template <
typename LO,
typename GO>
410 TEUCHOS_TEST_FOR_EXCEPTION(!buildConnectivityRun_,std::logic_error,
"DOFManager::getGIDFieldOffsets: cannot be called before "
411 "buildGlobalUnknowns has been called");
412 std::map<std::string,int>::const_iterator bitr = blockNameToID_.find(blockID);
413 if(bitr==blockNameToID_.end())
415 int bid=bitr->second;
416 if(fa_fps_[bid]!=Teuchos::null)
417 return fa_fps_[bid]->localOffsetsKokkos(fieldNum);
419 static const Kokkos::View<int*,PHX::Device> empty(
"panzer::DOFManager::getGIDFieldOffsetsKokkos() empty",0);
423 template <
typename LO,
typename GO>
426 gids = elementGIDs_[localElementID];
429 template <
typename LocalOrdinalT,
typename GlobalOrdinalT>
435 if(requireOrientations_){
441 std::vector<std::pair<FieldType,Teuchos::RCP<const FieldPattern>>> tmp;
442 for (std::size_t i=0; i < fieldPatterns_.size(); ++i)
443 tmp.push_back(std::make_pair(fieldTypes_[i],fieldPatterns_[i]));
447 connMngr_->buildConnectivity(*aggFieldPattern);
450 buildGlobalUnknowns(aggFieldPattern);
454 template <
typename LO,
typename GO>
459 typedef Tpetra::Map<LO, GO, Node> Map;
461 typedef Tpetra::Import<LO,GO,Node> Import;
464 typedef Tpetra::MultiVector<GO,LO,GO,Node> MultiVector;
466 PANZER_FUNC_TIME_MONITOR_DIFF(
"panzer::DOFManager::buildGlobalUnknowns",buildGlobalUnknowns);
473 "DOFManager::buildGlobalUnknowns: buildGlobalUnknowns cannot be called again "
474 "after buildGlobalUnknowns has been called");
478 if(getOrientationsRequired()) {
482 "DOFManager::buildGlobalUnknowns requires a geometric pattern including "
483 "the nodes when orientations are needed!");
489 ga_fp_ = geomPattern;
500 RCP<MultiVector> overlap_mv = Tpetra::createMultiVector<GO>(overlap_map,(size_t)numFields_);
503 auto non_overlap_pair = buildGlobalUnknowns_GUN(*tagged_overlap_mv,*overlap_mv);
514 fillGIDsFromOverlappedMV(ownedAccess,elementGIDs_,*overlap_map,*overlap_mv);
522 buildOverlapMapFromElements(neighborAccess);
525 Import imp_neighbor(non_overlap_map,overlap_map_neighbor);
528 Tpetra::createMultiVector<GO>(overlap_map_neighbor, (size_t)numFields_);
531 overlap_mv_neighbor->doImport(*non_overlap_mv, imp_neighbor,
534 fillGIDsFromOverlappedMV(neighborAccess, elementGIDs_,
535 *overlap_map_neighbor, *overlap_mv_neighbor);
544 panzer::Ordinal64 offset = 0xFFFFFFFFLL;
546 for (
size_t b = 0; b < blockOrder_.size(); ++b) {
547 const std::vector<LO> & myElements = connMngr_->getElementBlock(blockOrder_[b]);
548 for (std::size_t l = 0; l < myElements.size(); ++l) {
549 std::vector<GO> & localGIDs = elementGIDs_[myElements[l]];
550 for(std::size_t c=0;c<localGIDs.size();c++)
551 localGIDs[c] += offset;
556 for (
int j = 0; j < nvals.
size(); ++j)
563 PANZER_FUNC_TIME_MONITOR_DIFF(
"panzer::DOFManager::buildGlobalUnknowns::build_owned_vector",build_owned_vector);
565 typedef std::unordered_set<GO> HashTable;
566 HashTable isOwned, remainingOwned;
569 non_overlap_mv->sync_host();
570 tagged_non_overlap_mv->sync_host();
571 PHX::Device::fence();
572 auto nvals = non_overlap_mv->template getLocalView<PHX::Device>();
573 auto tagged_vals = tagged_non_overlap_mv->template getLocalView<PHX::Device>();
575 for (
size_t i = 0; i < nvals.extent(1); ++i) {
576 for (
size_t j = 0; j < nvals.extent(0); ++j) {
577 if (nvals(j,i) != -1) {
578 for(
GO offset=0;offset<tagged_vals(j,i);++offset)
579 isOwned.insert(nvals(j,i)+offset);
587 remainingOwned = isOwned;
590 for (
size_t b = 0; b < blockOrder_.size(); ++b) {
592 if(fa_fps_[b]==Teuchos::null)
595 const std::vector<LO> & myElements = connMngr_->getElementBlock(blockOrder_[b]);
597 for (
size_t l = 0; l < myElements.size(); ++l) {
598 const std::vector<GO> & localOrdering = elementGIDs_[myElements[l]];
601 for(std::size_t i=0;i<localOrdering.size();i++) {
603 if(isOwned.find(localOrdering[i])==isOwned.end())
607 std::pair<typename HashTable::iterator,bool> insertResult = hashTable.insert(localOrdering[i]);
608 if(insertResult.second) {
609 owned_.push_back(localOrdering[i]);
610 remainingOwned.erase(localOrdering[i]);
619 for(
typename HashTable::const_iterator itr=remainingOwned.begin();itr!=remainingOwned.end();itr++)
620 owned_.push_back(*itr);
622 if(owned_.size()!=isOwned.size()) {
623 out <<
"I'm about to hang because of unknown numbering failure ... sorry! (line = " << __LINE__ <<
")" << std::endl;
625 "DOFManager::buildGlobalUnkonwns: Failure because not all owned unknowns have been accounted for.");
636 PANZER_FUNC_TIME_MONITOR_DIFF(
"panzer::DOFManager::buildGlobalUnknowns::build_ghosted_array",BGA);
639 typedef std::unordered_set<GO> HashTable;
641 for (std::size_t i = 0; i < owned_.size(); i++)
642 hashTable.insert(owned_[i]);
647 std::vector<ElementBlockAccess> blockAccessVec;
651 for (std::size_t a = 0; a < blockAccessVec.size(); ++a)
655 for (
size_t b = 0; b < blockOrder_.size(); ++b)
657 if (fa_fps_[b] == Teuchos::null)
659 const std::vector<LO>& myElements =
661 for (
size_t l = 0; l < myElements.size(); ++l)
663 const std::vector<GO>& localOrdering = elementGIDs_[myElements[l]];
666 for (std::size_t i = 0; i < localOrdering.size(); ++i)
668 std::pair<typename HashTable::iterator, bool> insertResult =
669 hashTable.insert(localOrdering[i]);
673 if(insertResult.second)
674 ghosted_.push_back(localOrdering[i]);
681 buildConnectivityRun_ =
true;
684 if(requireOrientations_) {
685 PANZER_FUNC_TIME_MONITOR_DIFF(
"panzer::DOFManager::buildGlobalUnknowns::build_orientation",BO);
686 buildUnknownsOrientation();
691 PANZER_FUNC_TIME_MONITOR_DIFF(
"panzer::DOFManager::buildGlobalUnknowns::build_local_ids_from_owned_and_ghosted",BLOFOG);
692 this->buildLocalIdsFromOwnedAndGhostedElements();
695 PANZER_FUNC_TIME_MONITOR_DIFF(
"panzer::DOFManager::buildGlobalUnknowns::build_local_ids",BLI);
696 this->buildLocalIds();
700 template <
typename LO,
typename GO>
701 std::pair<Teuchos::RCP<Tpetra::MultiVector<GO,LO,GO,panzer::TpetraNodeType> >,
704 Tpetra::MultiVector<GO,LO,GO,panzer::TpetraNodeType> & overlap_mv)
const
708 typedef Tpetra::Map<LO, GO, Node> Map;
710 typedef Tpetra::Export<LO,GO,Node> Export;
711 typedef Tpetra::Import<LO,GO,Node> Import;
714 typedef Tpetra::MultiVector<GO,LO,GO,Node> MultiVector;
716 PANZER_FUNC_TIME_MONITOR_DIFF(
"panzer::DOFManager::buildGlobalUnknowns_GUN",BGU_GUN);
728 PANZER_FUNC_TIME_MONITOR_DIFF(
"panzer::DOFManager::buildGlobalUnknowns_GUN::line_04 createOneToOne",GUN04);
730 GreedyTieBreak<LO,GO> greedy_tie_break;
731 non_overlap_map = Tpetra::createOneToOne<LO,GO,Node>(overlap_map, greedy_tie_break);
736 HashTieBreak<LO,GO> tie_break;
737 non_overlap_map = Tpetra::createOneToOne<LO,GO,Node>(overlap_map,tie_break);
747 PANZER_FUNC_TIME_MONITOR_DIFF(
"panzer::DOFManager::buildGlobalUnknowns_GUN::line_05 alloc_unique_mv",GUN05);
749 tagged_non_overlap_mv = Tpetra::createMultiVector<GO>(non_overlap_map,(size_t)numFields_);
760 PANZER_FUNC_TIME_MONITOR_DIFF(
"panzer::DOFManager::buildGlobalUnknowns_GUN::line_06 export",GUN06);
762 exp =
rcp(
new Export(overlap_map,non_overlap_map));
764 #ifdef PANZER_DOFMGR_REQUIRE_CUDA
768 imp =
rcp(
new Import(non_overlap_map,overlap_map));
773 tagged_non_overlap_mv->doExport(tagged_overlap_mv,*exp,Tpetra::ABSMAX);
777 non_overlap_mv =
rcp(
new MultiVector(*tagged_non_overlap_mv,
Teuchos::Copy));
787 PANZER_FUNC_TIME_MONITOR_DIFF(
"panzer::DOFManager::buildGlobalUnknowns_GUN::line_07-09 local_count",GUN07_09);
788 auto values = non_overlap_mv->getLocalViewDevice();
789 auto mv_size = values.extent(0);
791 PHX::Device::fence();
801 PANZER_FUNC_TIME_MONITOR_DIFF(
"panzer::DOFManager::buildGlobalUnknowns_GUN::line_10 prefix_sum",GUN_10);
805 Teuchos::scan<int, GO> (*getComm(),
Teuchos::REDUCE_SUM,
static_cast<GO> (localsum), Teuchos::outArg (scanResult));
806 myOffset = scanResult - localsum;
819 PANZER_FUNC_TIME_MONITOR_DIFF(
"panzer::DOFManager::buildGlobalUnknowns_GUN::line_13-21 gid_assignment",GUN13_21);
821 auto editnonoverlap = non_overlap_mv->getLocalViewHost();
822 for(
size_t i=0; i<non_overlap_mv->getLocalLength(); ++i){
823 for(
int j=0; j<numFields_; ++j){
824 if(editnonoverlap(i,j)!=0){
826 int ndof = Teuchos::as<int>(editnonoverlap(i,j));
827 editnonoverlap(i,j)=myOffset+which_id;
831 editnonoverlap(i,j)=-1;
836 non_overlap_mv->modify_host();
837 non_overlap_mv->sync_device();
838 PHX::Device::fence();
850 PANZER_FUNC_TIME_MONITOR_DIFF(
"panzer::DOFManager::buildGlobalUnknowns_GUN::line_23 final_import",GUN23);
852 #ifdef PANZER_DOFMGR_REQUIRE_CUDA
853 overlap_mv.doImport(*non_overlap_mv,*imp,Tpetra::REPLACE);
856 overlap_mv.doImport(*non_overlap_mv,*exp,Tpetra::REPLACE);
863 return std::make_pair(non_overlap_mv,tagged_non_overlap_mv);
866 template <
typename LO,
typename GO>
872 typedef Tpetra::Map<LO, GO, Node> Map;
873 typedef Tpetra::MultiVector<GO,LO,GO,Node> MultiVector;
875 PANZER_FUNC_TIME_MONITOR_DIFF(
"panzer::DOFManager::buildTaggedMultiVector",BTMV);
880 for (
size_t b = 0; b < blockOrder_.size(); ++b) {
881 std::vector<std::tuple< int, panzer::FieldType, RCP<const panzer::FieldPattern> > > faConstruct;
885 for (
size_t i = 0; i < fieldAIDOrder_.size(); ++i) {
886 int looking = fieldAIDOrder_[i];
889 std::vector<int>::const_iterator reu = std::find(blockToAssociatedFP_[b].begin(), blockToAssociatedFP_[b].
end(), looking);
890 if(!(reu==blockToAssociatedFP_[b].
end())){
891 faConstruct.push_back(std::make_tuple(i, fieldTypes_[fieldAIDOrder_[i]], fieldPatterns_[fieldAIDOrder_[i]]));
896 if(faConstruct.size()>0) {
900 int gidsInBlock = fa_fps_[fa_fps_.size()-1]->numberIds();
901 elementBlockGIDCount_.push_back(gidsInBlock);
904 fa_fps_.push_back(Teuchos::null);
905 elementBlockGIDCount_.push_back(0);
909 RCP<const Map> overlapmap = buildOverlapMapFromElements(ownedAccess);
910 PHX::Device::fence();
917 PANZER_FUNC_TIME_MONITOR_DIFF(
"panzer::DOFManager::buildTaggedMultiVector::allocate_tagged_multivector",ATMV);
919 overlap_mv = Tpetra::createMultiVector<GO>(overlapmap,(size_t)numFields_);
920 overlap_mv->putScalar(0);
921 PHX::Device::fence();
929 PANZER_FUNC_TIME_MONITOR_DIFF(
"panzer::DOFManager::buildTaggedMultiVector::fill_tagged_multivector",FTMV);
932 std::vector<int> working(overlap_mv->getNumVectors());
933 overlap_mv->sync_host();
934 PHX::Device::fence();
935 auto edittwoview_host = overlap_mv->getLocalViewHost();
936 for (
size_t b = 0; b < blockOrder_.size(); ++b) {
938 if(fa_fps_[b]==Teuchos::null)
941 const std::vector<LO> &
numFields= fa_fps_[b]->numFieldsPerId();
942 const std::vector<LO> & fieldIds= fa_fps_[b]->fieldIds();
943 const std::vector<LO> & myElements = connMngr_->getElementBlock(blockOrder_[b]);
944 for (
size_t l = 0; l < myElements.size(); ++l) {
945 LO connSize = connMngr_->getConnectivitySize(myElements[l]);
946 const GO * elmtConn = connMngr_->getConnectivity(myElements[l]);
948 for (
int c = 0; c < connSize; ++c) {
949 size_t lid = overlapmap->getLocalElement(elmtConn[c]);
951 for(std::size_t i=0;i<working.size();i++)
953 for (
int n = 0; n < numFields[c]; ++n) {
954 int whichField = fieldIds[offset];
957 working[whichField]++;
960 for(std::size_t i=0;i<working.size();i++) {
961 auto current = edittwoview_host(lid,i);
962 edittwoview_host(lid,i) = (current > working[i]) ? current : working[i];
968 overlap_mv->modify_host();
969 overlap_mv->sync_device();
970 PHX::Device::fence();
983 template <
typename LO,
typename GO>
988 while(!found && (
size_t)ind<fieldStringOrder_.size()){
989 if(fieldStringOrder_[ind]==
string)
1002 template <
typename LO,
typename GO>
1005 fieldOrder.resize(fieldStringOrder_.size());
1006 for (
size_t i = 0; i < fieldStringOrder_.size(); ++i)
1007 fieldOrder[i]=fieldStringOrder_[i];
1010 template <
typename LO,
typename GO>
1013 std::map<std::string,int>::const_iterator fitr = fieldNameToAID_.find(field);
1014 if(fitr==fieldNameToAID_.end()) {
1015 std::stringstream ss;
1016 printFieldInformation(ss);
1017 TEUCHOS_TEST_FOR_EXCEPTION(
true,std::logic_error,
"DOFManager::fieldInBlock: invalid field name. DOF information is:\n"+ss.str());
1019 std::map<std::string,int>::const_iterator bitr = blockNameToID_.find(block);
1020 if(bitr==blockNameToID_.end()) {
1021 std::stringstream ss;
1022 printFieldInformation(ss);
1023 TEUCHOS_TEST_FOR_EXCEPTION(
true,std::logic_error,
"DOFManager::fieldInBlock: invalid block name. DOF information is:\n"+ss.str());
1025 int fid=fitr->second;
1026 int bid=bitr->second;
1029 for (
size_t i = 0; i < blockToAssociatedFP_[bid].size(); ++i) {
1030 if(blockToAssociatedFP_[bid][i]==fid){
1039 template <
typename LO,
typename GO>
1042 TEUCHOS_TEST_FOR_EXCEPTION(!buildConnectivityRun_,std::logic_error,
"DOFManager::getBlockFieldNumbers: BuildConnectivity must be run first.");
1044 std::map<std::string,int>::const_iterator bitr = blockNameToID_.find(blockId);
1045 if(bitr==blockNameToID_.end())
1047 int bid=bitr->second;
1050 if(fa_fps_[bid]!=Teuchos::null)
1051 return fa_fps_[bid]->fieldIds();
1054 static std::vector<int> empty;
1058 template <
typename LO,
typename GO>
1059 const std::pair<std::vector<int>,std::vector<int> > &
1062 TEUCHOS_TEST_FOR_EXCEPTION(!buildConnectivityRun_,std::logic_error,
"DOFManager::getGIDFieldOffsets_closure: BuildConnectivity must be run first.");
1063 std::map<std::string,int>::const_iterator bitr = blockNameToID_.find(blockId);
1064 if(bitr==blockNameToID_.end())
1068 if(fa_fps_[bitr->second]!=Teuchos::null)
1069 return fa_fps_[bitr->second]->localOffsets_closure(fieldNum, subcellDim, subcellId);
1071 static std::pair<std::vector<int>,std::vector<int> > empty;
1075 template <
typename LO,
typename GO>
1079 if(indices.size()!=isOwned.size())
1080 isOwned.resize(indices.size(),
false);
1081 typename std::vector<GO>::const_iterator endOf = owned_.end();
1082 for (std::size_t i = 0; i < indices.size(); ++i) {
1083 isOwned[i] = ( std::find(owned_.begin(), owned_.end(), indices[i])!=endOf );
1087 template <
typename LO,
typename GO>
1091 "DOFManager::setFieldOrder: setFieldOrder cannot be called after "
1092 "buildGlobalUnknowns has been called");
1093 if(validFieldOrder(fieldOrder)){
1094 fieldStringOrder_=fieldOrder;
1096 for (
size_t i = 0; i < fieldStringOrder_.size(); ++i) {
1097 fieldAIDOrder_[i]=fieldNameToAID_[fieldStringOrder_[i]];
1106 template <
typename LO,
typename GO>
1109 if(fieldStringOrder_.size()!=proposed_fieldOrder.size())
1113 for (
size_t i = 0; i < fieldStringOrder_.size(); ++i) {
1115 for (
size_t j = 0; j < proposed_fieldOrder.size(); ++j) {
1116 if(fieldStringOrder_[i]==proposed_fieldOrder[j])
1125 template <
typename LO,
typename GO>
1129 if(num>=(
int)fieldStringOrder_.size())
1131 return fieldStringOrder_[num];
1139 template <
typename LO,
typename GO>
1142 orientation_.clear();
1144 std::vector<std::string> elementBlockIds;
1145 connMngr_->getElementBlockIds(elementBlockIds);
1148 std::size_t myElementCount = 0;
1149 for(std::vector<std::string>::const_iterator blockItr=elementBlockIds.begin(); blockItr!=elementBlockIds.end();++blockItr)
1150 myElementCount += connMngr_->getElementBlock(*blockItr).size();
1153 orientation_.resize(myElementCount);
1156 for(std::vector<std::string>::const_iterator blockItr=elementBlockIds.begin();
1157 blockItr!=elementBlockIds.end();++blockItr) {
1158 const std::string & blockName = *blockItr;
1161 std::map<std::string,int>::const_iterator fap = blockNameToID_.find(blockName);
1162 if(fap==blockNameToID_.end()) {
1166 int bid=fap->second;
1168 if(fa_fps_[bid]==Teuchos::null)
1175 std::vector<std::pair<int,int> > topEdgeIndices;
1179 std::vector<std::vector<int> > topFaceIndices;
1180 if(ga_fp_->getDimension()==3)
1184 const std::vector<LO> & elmts = connMngr_->getElementBlock(blockName);
1185 for(std::size_t e=0;e<elmts.size();e++) {
1187 std::vector<signed char> & eOrientation = orientation_[elmts[e]];
1192 eOrientation.resize(fieldPattern.
numberIds());
1193 for(std::size_t s=0;s<eOrientation.size();s++)
1194 eOrientation[s] = 1;
1197 LO connSz = connMngr_->getConnectivitySize(elmts[e]);
1198 const GO * connPtr = connMngr_->getConnectivity(elmts[e]);
1199 const std::vector<GO> connectivity(connPtr,connPtr+connSz);
1204 if(ga_fp_->getDimension()==3)
1210 template <
typename LO,
typename GO>
1214 "DOFManager::getElementOrientations: Orientations were not constructed!");
1216 const std::vector<signed char> & local_o = orientation_[localElmtId];
1217 gidsOrientation.resize(local_o.size());
1218 for(std::size_t i=0;i<local_o.size();i++) {
1219 gidsOrientation[i] = double(local_o[i]);
1223 template <
typename LocalOrdinalT,
typename GlobalOrdinalT>
1228 connMngr_ = Teuchos::null;
1231 orientation_.clear();
1233 elementGIDs_.clear();
1236 elementBlockGIDCount_.clear();
1241 template <
typename LocalOrdinalT,
typename GlobalOrdinalT>
1244 std::map<std::string,int>::const_iterator bitr = blockNameToID_.find(blockId);
1245 if(bitr==blockNameToID_.end())
1247 return bitr->second;
1250 template <
typename LocalOrdinalT,
typename GlobalOrdinalT>
1253 os <<
"DOFManager Field Information: " << std::endl;
1258 for(std::size_t i=0;i<blockOrder_.size();i++) {
1259 os <<
" Element Block = " << blockOrder_[i] << std::endl;
1262 const std::vector<int> & fieldIds = blockToAssociatedFP_[i];
1263 for(std::size_t f=0;f<fieldIds.size();f++)
1264 os <<
" \"" << getFieldString(fieldIds[f]) <<
"\" is field ID " << fieldIds[f] << std::endl;
1268 template <
typename LO,
typename GO>
1273 PANZER_FUNC_TIME_MONITOR(
"panzer::DOFManager::builderOverlapMapFromElements");
1280 std::set<GO> overlapset;
1281 for (
size_t i = 0; i < blockOrder_.size(); ++i) {
1282 const std::vector<LO> & myElements = access.
getElementBlock(blockOrder_[i]);
1283 for (
size_t e = 0; e < myElements.size(); ++e) {
1284 LO connSize = connMngr_->getConnectivitySize(myElements[e]);
1285 const GO * elmtConn = connMngr_->getConnectivity(myElements[e]);
1286 for (
int k = 0; k < connSize; ++k) {
1287 overlapset.insert(elmtConn[k]);
1293 for (
typename std::set<GO>::const_iterator itr = overlapset.begin(); itr!=overlapset.end(); ++itr) {
1299 return Tpetra::createNonContigMap<LO,GO>(overlapVector,getComm());
1302 template <
typename LO,
typename GO>
1305 std::vector<std::vector< GO > > & elementGIDs,
1306 const Tpetra::Map<LO,GO,panzer::TpetraNodeType> & overlapmap,
1307 const Tpetra::MultiVector<GO,LO,GO,panzer::TpetraNodeType> & const_overlap_mv)
const
1312 auto overlap_mv =
const_cast<Tpetra::MultiVector<GO,LO,GO,panzer::TpetraNodeType>&
>(const_overlap_mv);
1313 overlap_mv.sync_host();
1314 PHX::Device::fence();
1315 const auto twoview_host = overlap_mv.getLocalViewHost();
1319 for (
size_t b = 0; b < blockOrder_.size(); ++b) {
1320 const std::vector<LO> & myElements = access.
getElementBlock(blockOrder_[b]);
1322 if(fa_fps_[b]==Teuchos::null) {
1324 for (
size_t l = 0; l < myElements.size(); ++l) {
1325 LO thisID=myElements[l];
1326 if(elementGIDs.size()<=(size_t)thisID)
1327 elementGIDs.resize(thisID+1);
1332 const std::vector<int> &
numFields= fa_fps_[b]->numFieldsPerId();
1333 const std::vector<int> & fieldIds= fa_fps_[b]->fieldIds();
1336 for (
size_t l = 0; l < myElements.size(); ++l) {
1337 LO connSize = connMngr_->getConnectivitySize(myElements[l]);
1338 const GO * elmtConn = connMngr_->getConnectivity(myElements[l]);
1339 std::vector<GO> localOrdering;
1341 for (
int c = 0; c < connSize; ++c) {
1342 size_t lid = overlapmap.getLocalElement(elmtConn[c]);
1343 std::vector<int> dofsPerField(numFields_,0);
1344 for (
int n = 0; n < numFields[c]; ++n) {
1345 int whichField = fieldIds[offset];
1349 localOrdering.push_back(twoview_host(lid,whichField)+dofsPerField[whichField]);
1351 dofsPerField[whichField]++;
1354 LO thisID=myElements[l];
1355 if(elementGIDs.size()<=(size_t)thisID){
1356 elementGIDs.resize(thisID+1);
1358 elementGIDs[thisID]=localOrdering;
1363 template <
typename LO,
typename GO>
1366 std::vector<std::vector<LO> > elementLIDs(elementGIDs_.size());
1368 std::vector<GO> ownedAndGhosted;
1369 this->getOwnedAndGhostedIndices(ownedAndGhosted);
1372 std::unordered_map<GO,LO> hashMap;
1373 for(std::size_t i = 0; i < ownedAndGhosted.size(); ++i)
1374 hashMap[ownedAndGhosted[i]] = i;
1376 for (std::size_t i = 0; i < elementGIDs_.size(); ++i) {
1377 const std::vector<GO>& gids = elementGIDs_[i];
1378 std::vector<LO>&
lids = elementLIDs[i];
1379 lids.resize(gids.size());
1380 for (std::size_t g = 0; g < gids.size(); ++g)
1381 lids[g] = hashMap[gids[g]];
1384 this->setLocalIds(elementLIDs);
int getNumOwned() const
Get the number of indices owned by this processor.
void buildLocalIdsFromOwnedAndGhostedElements()
void computeCellEdgeOrientations(const std::vector< std::pair< int, int > > &topEdgeIndices, const std::vector< GlobalOrdinalT > &topology, const FieldPattern &fieldPattern, std::vector< signed char > &orientation)
Sums all entries of a Rank 2 Kokkos View.
void getFieldOrder(std::vector< std::string > &fieldOrder) const
const std::vector< int > & getBlockFieldNumbers(const std::string &blockId) const
basic_FancyOStream & setShowProcRank(const bool showProcRank)
const std::vector< int > & getGIDFieldOffsets(const std::string &blockID, int fieldNum) const
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
void setConnManager(const Teuchos::RCP< ConnManager< LO, GO > > &connMngr, MPI_Comm mpiComm)
Adds a Connection Manager that will be associated with this DOFManager.
void printFieldInformation(std::ostream &os) const
const Kokkos::View< const int *, PHX::Device > getGIDFieldOffsetsKokkos(const std::string &blockID, int fieldNum) const
void getElementOrientation(LocalOrdinalT localElmtId, std::vector< double > &gidsOrientation) const
Get a vector containg the orientation of the GIDs relative to the neighbors.
FieldType
The type of discretization to use for a field pattern.
int getNumGhosted() const
Get the number of indices ghosted for this processor.
void ownedIndices(const std::vector< GlobalOrdinalT > &indices, std::vector< bool > &isOwned) const
virtual int numberIds() const
Teuchos::RCP< ConnManager< LocalOrdinalT, GlobalOrdinalT > > resetIndices()
Reset the indices for this DOF manager.
int getNumOwnedAndGhosted() const
Get the number of owned and ghosted indices for this processor.
Kokkos::View< const LO **, PHX::Device > lids
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Teuchos::RCP< const FieldPattern > getFieldPattern(const std::string &name) const
Find a field pattern stored for a particular block and field number. This will retrive the pattern ad...
void getElementGIDs(LO localElementID, std::vector< GO > &gids, const std::string &blockIdHint="") const
get associated GIDs for a given local element
int addField(const std::string &str, const Teuchos::RCP< const FieldPattern > &pattern, const panzer::FieldType &type=panzer::FieldType::CG)
Add a field to the DOF manager.
void getGhostedIndices(std::vector< GlobalOrdinalT > &indices) const
Get the set of indices ghosted for this processor.
std::pair< Teuchos::RCP< Tpetra::MultiVector< GO, LO, GO, panzer::TpetraNodeType > >, Teuchos::RCP< Tpetra::MultiVector< GO, LO, GO, panzer::TpetraNodeType > > > buildGlobalUnknowns_GUN(const Tpetra::MultiVector< GO, LO, GO, panzer::TpetraNodeType > &tagged_overlap_mv, Tpetra::MultiVector< GO, LO, GO, panzer::TpetraNodeType > &overlap_mv) const
void setFieldOrder(const std::vector< std::string > &fieldOrder)
Kokkos::Compat::KokkosDeviceWrapperNode< PHX::Device > TpetraNodeType
void buildUnknownsOrientation()
basic_FancyOStream & setOutputToRootOnly(const int rootRank)
void buildGlobalUnknowns()
builds the global unknowns array
void computePatternFaceIndices(const FieldPattern &pattern, std::vector< std::vector< int > > &faceIndices)
bool validFieldOrder(const std::vector< std::string > &proposed_fieldOrder)
Teuchos::RCP< const Tpetra::Map< LO, GO, panzer::TpetraNodeType > > buildOverlapMapFromElements(const ElementBlockAccess &access) const
void getOwnedAndGhostedIndices(std::vector< GlobalOrdinalT > &indices) const
Get the set of owned and ghosted indices for this processor.
const std::string & getFieldString(int num) const
Reverse lookup of the field string from a field number.
void computeCellFaceOrientations(const std::vector< std::pair< int, int > > &topEdgeIndices, const std::vector< GlobalOrdinalT > &topology, const FieldPattern &fieldPattern, std::vector< signed char > &orientation)
Teuchos::RCP< Tpetra::MultiVector< GO, LO, GO, panzer::TpetraNodeType > > buildTaggedMultiVector(const ElementBlockAccess &access)
void push_back(const value_type &x)
PHX::MDField< ScalarT, panzer::Cell, panzer::BASIS > field
A field to which we'll contribute, or in which we'll store, the result of computing this integral...
bool fieldInBlock(const std::string &field, const std::string &block) const
void getOwnedIndices(std::vector< GlobalOrdinalT > &indices) const
Get the set of indices owned by this processor.
int getNumFields() const
gets the number of fields
std::size_t blockIdToIndex(const std::string &blockId) const
const std::pair< std::vector< int >, std::vector< int > > & getGIDFieldOffsets_closure(const std::string &blockId, int fieldNum, int subcellDim, int subcellId) const
Use the field pattern so that you can find a particular field in the GIDs array. This version lets yo...
int getFieldNum(const std::string &string) const
Get the number used for access to this field.
virtual const std::vector< int > & getSubcellIndices(int dim, int cellIndex) const =0
#define TEUCHOS_ASSERT(assertion_test)
void fillGIDsFromOverlappedMV(const ElementBlockAccess &access, std::vector< std::vector< GO > > &elementGIDs, const Tpetra::Map< LO, GO, panzer::TpetraNodeType > &overlapmap, const Tpetra::MultiVector< GO, LO, GO, panzer::TpetraNodeType > &overlap_mv) const
const std::vector< LO > & getElementBlock(const std::string &eBlock) const
void computePatternEdgeIndices(const FieldPattern &pattern, std::vector< std::pair< int, int > > &edgeIndices)