43 #ifndef __Panzer_Filtered_UniqueGlobalIndexer_impl_hpp__ 44 #define __Panzer_Filtered_UniqueGlobalIndexer_impl_hpp__ 46 #include <unordered_set> 48 #include "PanzerDofMgr_config.hpp" 51 #include "Tpetra_Map.hpp" 52 #include "Tpetra_Import.hpp" 53 #include "Tpetra_Export.hpp" 54 #include "Tpetra_Vector.hpp" 58 template <
typename LocalOrdinalT,
typename GlobalOrdinalT>
68 template<
typename LocalOrdinalT,
typename GlobalOrdinalT>
73 const std::vector<GlobalOrdinalT>& filtered)
75 typedef GlobalOrdinalT GO;
76 typedef LocalOrdinalT LO;
78 typedef Tpetra::Map<LO, GO, Node> Map;
79 typedef Tpetra::Vector<GO,LO,GO,Node> Vector;
80 typedef Tpetra::Export<LO,GO,Node> Export;
84 using HashTable = std::unordered_set<GlobalOrdinalT>;
92 vector<GlobalOrdinalT> baseOwned, baseGhosted;
93 base_->getOwnedIndices(baseOwned);
94 base_->getGhostedIndices(baseGhosted);
97 = Tpetra::createNonContigMap<LO,GO>(baseOwned,getComm());
99 = Tpetra::createNonContigMap<LO,GO>(baseGhosted,getComm());
101 Vector ownedFiltered(ownedMap);
102 Vector ghostedFiltered(ghostedMap);
104 ownedFiltered.putScalar(0.0);
105 ghostedFiltered.putScalar(0.0);
107 for(GlobalOrdinalT f : filtered) {
108 bool isOwned = std::find(baseOwned.begin(),baseOwned.end(),f)!=baseOwned.end();
109 bool isGhosted = std::find(baseGhosted.begin(),baseGhosted.end(),f)!=baseGhosted.end();
112 ownedFiltered.replaceGlobalValue(f,1.0);
114 ghostedFiltered.replaceGlobalValue(f,1.0);
118 Export exporter(ghostedMap,ownedMap);
119 ownedFiltered.doExport(ghostedFiltered, exporter, Tpetra::ADD);
124 HashTable filteredHash;
125 for(
int i(0); i < data.
size(); ++i) {
127 filteredHash.insert(baseOwned[i]);
135 for (
size_t i(0); i < baseOwned.size(); ++i)
137 auto itr = filteredHash.find(baseOwned[i]);
138 if (itr == filteredHash.end())
139 owned_.push_back(baseOwned[i]);
141 ghosted_.push_back(baseOwned[i]);
143 ghosted_.insert(ghosted_.end(), baseGhosted.begin(), baseGhosted.end());
147 this->buildLocalIds();
150 template <
typename LocalOrdinalT,
typename GlobalOrdinalT>
157 typedef GlobalOrdinalT GO;
158 typedef LocalOrdinalT LO;
160 typedef Tpetra::Map<LO, GO, Node> Map;
161 typedef Tpetra::Vector<GO,LO,GO,Node> Vector;
162 typedef Tpetra::Import<LO,GO,Node> Import;
164 std::vector<GlobalOrdinalT> ownedIndices;
165 std::vector<GlobalOrdinalT> ghostedIndices;
168 getOwnedIndices(ownedIndices);
169 getOwnedAndGhostedIndices(ghostedIndices);
172 = Tpetra::createNonContigMap<LO,GO>(ownedIndices,getComm());
174 = Tpetra::createNonContigMap<LO,GO>(ghostedIndices,getComm());
178 Vector ownedActive(ownedMap);
179 ownedActive.putScalar(1);
182 Vector ghostedActive(ghostedMap);
183 ghostedActive.putScalar(0);
187 Import importer(ownedMap,ghostedMap);
188 ghostedActive.doImport(ownedActive,importer,Tpetra::INSERT);
194 indicator.insert(indicator.end(),data.
begin(),data.
end());
197 template <
typename LocalOrdinalT,
typename GlobalOrdinalT>
205 std::vector<int> indicators;
206 getOwnedAndGhostedNotFilteredIndicator(indicators);
209 std::vector<GlobalOrdinalT> ghostedIndices;
210 getOwnedAndGhostedIndices(ghostedIndices);
213 for(std::size_t i=0;i<indicators.size();i++) {
215 indices.push_back(ghostedIndices[i]);
219 template <
typename LocalOrdinalT,
typename GlobalOrdinalT>
222 ownedIndices(
const std::vector<GlobalOrdinalT> & indices,std::vector<bool> & isOwned)
const 225 if(indices.size()!=isOwned.size())
226 isOwned.resize(indices.size(),
false);
227 typename std::vector<GlobalOrdinalT>::const_iterator endOf = owned_.end();
228 for (std::size_t i = 0; i < indices.size(); ++i) {
229 isOwned[i] = ( std::find(owned_.begin(), owned_.end(), indices[i])!=endOf );
void initialize(const Teuchos::RCP< const UniqueGlobalIndexer< LocalOrdinalT, GlobalOrdinalT > > &ugi, const std::vector< GlobalOrdinalT > &filteredIndices)
Filtered_UniqueGlobalIndexer()
void getOwnedAndGhostedNotFilteredIndicator(std::vector< int > &indicator) const
Kokkos::Compat::KokkosDeviceWrapperNode< PHX::Device > TpetraNodeType
virtual void ownedIndices(const std::vector< GlobalOrdinalT > &indices, std::vector< bool > &isOwned) const
void getFilteredOwnedAndGhostedIndices(std::vector< GlobalOrdinalT > &indices) const