50 #include "Kokkos_DynRankView.hpp"
57 #include "Teuchos_FancyOStream.hpp"
59 namespace panzer_stk {
71 bool operator() (stk::mesh::Entity a,stk::mesh::Entity b)
72 {
return mesh_->elementLocalId(a) <
mesh_->elementLocalId(b);}
78 template <
typename GO>
80 : stkMeshDB_(stkMeshDB), ownedElementCount_(0)
84 template <
typename GO>
92 template <
typename GO>
95 elements_ = Teuchos::null;
97 elementBlocks_.clear();
98 elmtLidToConn_.clear();
100 elmtToAssociatedElmts_.clear();
103 template <
typename GO>
106 clearLocalElementMapping();
110 elements_ =
Teuchos::rcp(
new std::vector<stk::mesh::Entity>);
113 std::vector<std::string> blockIds;
114 stkMeshDB_->getElementBlockNames(blockIds);
116 std::size_t blockIndex=0;
117 for(std::vector<std::string>::const_iterator idItr=blockIds.begin();
118 idItr!=blockIds.end();++idItr,++blockIndex) {
119 std::string blockId = *idItr;
122 std::vector<stk::mesh::Entity> blockElmts;
123 stkMeshDB_->getMyElements(blockId,blockElmts);
126 elements_->insert(elements_->end(),blockElmts.begin(),blockElmts.end());
129 elementBlocks_[blockId] =
Teuchos::rcp(
new std::vector<LocalOrdinal>);
130 for(std::size_t i=0;i<blockElmts.size();i++)
131 elementBlocks_[blockId]->push_back(stkMeshDB_->elementLocalId(blockElmts[i]));
134 ownedElementCount_ = elements_->size();
137 for(std::vector<std::string>::const_iterator idItr=blockIds.begin();
138 idItr!=blockIds.end();++idItr,++blockIndex) {
139 std::string blockId = *idItr;
142 std::vector<stk::mesh::Entity> blockElmts;
143 stkMeshDB_->getNeighborElements(blockId,blockElmts);
146 elements_->insert(elements_->end(),blockElmts.begin(),blockElmts.end());
149 neighborElementBlocks_[blockId] =
Teuchos::rcp(
new std::vector<LocalOrdinal>);
150 for(std::size_t i=0;i<blockElmts.size();i++)
151 neighborElementBlocks_[blockId]->push_back(stkMeshDB_->elementLocalId(blockElmts[i]));
155 std::sort(elements_->begin(),elements_->end(),
LocalIdCompare(stkMeshDB_));
159 elmtLidToConn_.clear();
160 elmtLidToConn_.resize(elements_->size(),0);
163 connSize_.resize(elements_->size(),0);
166 template <
typename GO>
174 GlobalOrdinal maxNodeId = stkMeshDB_->getMaxEntityId(stkMeshDB_->getNodeRank());
175 GlobalOrdinal maxEdgeId = stkMeshDB_->getMaxEntityId(stkMeshDB_->getEdgeRank());
176 GlobalOrdinal maxFaceId = stkMeshDB_->getMaxEntityId(stkMeshDB_->getFaceRank());
198 edgeOffset = nodeOffset+(maxNodeId+1)*nodeIdCnt;
199 faceOffset = edgeOffset+(maxEdgeId+1)*edgeIdCnt;
200 cellOffset = faceOffset+(maxFaceId+1)*faceIdCnt;
204 && edgeOffset <= faceOffset
205 && faceOffset <= cellOffset);
208 template <
typename GO>
217 stk::mesh::BulkData& bulkData = *stkMeshDB_->getBulkData();
218 const stk::mesh::EntityRank rank =
static_cast<stk::mesh::EntityRank
>(subcellRank);
219 const size_t num_rels = bulkData.num_connectivity(element, rank);
220 stk::mesh::Entity
const* relations = bulkData.begin(element, rank);
221 for(std::size_t sc=0; sc<num_rels; ++sc) {
222 stk::mesh::Entity subcell = relations[sc];
226 connectivity_.push_back(offset+idCnt*(bulkData.identifier(subcell)-1)+i);
233 template <
typename GO>
238 LocalOrdinal elmtLID = stkMeshDB_->elementLocalId(element);
240 const std::vector<int> & subCellIndices = fp.
getSubcellIndices(subcellRank,subcellId);
243 for(std::size_t i=0;i<subCellIndices.size();i++) {
244 conn[subCellIndices[i]] = offset+subCellIndices.size()*(newId-1)+i;
248 template <
typename GO>
251 #ifdef HAVE_EXTRA_TIMERS
253 RCP<Teuchos::TimeMonitor> tM =
rcp(
new TimeMonitor(*TimeMonitor::getNewTimer(std::string(
"panzer_stk::STKConnManager::buildConnectivity"))));
256 stk::mesh::BulkData& bulkData = *stkMeshDB_->getBulkData();
260 buildLocalElementMapping();
266 LocalOrdinal nodeIdCnt=0, edgeIdCnt=0, faceIdCnt=0, cellIdCnt=0;
267 GlobalOrdinal nodeOffset=0, edgeOffset=0, faceOffset=0, cellOffset=0;
268 buildOffsetsAndIdCounts(fp, nodeIdCnt, edgeIdCnt, faceIdCnt, cellIdCnt,
269 nodeOffset, edgeOffset, faceOffset, cellOffset);
277 for(std::size_t elmtLid=0;elmtLid!=elements_->size();++elmtLid) {
279 stk::mesh::Entity element = (*elements_)[elmtLid];
282 elmtLidToConn_[elmtLid] = connectivity_.size();
285 numIds += addSubcellConnectivities(element,stkMeshDB_->getNodeRank(),nodeIdCnt,nodeOffset);
286 numIds += addSubcellConnectivities(element,stkMeshDB_->getEdgeRank(),edgeIdCnt,edgeOffset);
287 numIds += addSubcellConnectivities(element,stkMeshDB_->getFaceRank(),faceIdCnt,faceOffset);
293 connectivity_.push_back(cellOffset+cellIdCnt*(bulkData.identifier(element)-1));
298 connSize_[elmtLid] = numIds;
301 applyPeriodicBCs( fp, nodeOffset, edgeOffset, faceOffset, cellOffset);
306 if (hasAssociatedNeighbors())
307 applyInterfaceConditions();
310 template <
typename GO>
314 stk::mesh::Entity element = (*elements_)[localElmtId];
316 return stkMeshDB_->containingBlockId(element);
319 template <
typename GO>
326 #ifdef HAVE_EXTRA_TIMERS
328 RCP<Teuchos::TimeMonitor> tM =
rcp(
new TimeMonitor(*TimeMonitor::getNewTimer(std::string(
"panzer_stk::STKConnManager::applyPeriodicBCs"))));
332 = stkMeshDB_->getPeriodicNodePairing();
335 = matchedValues.first;
337 = matchedValues.second;
340 if(matchedNodes==Teuchos::null)
return;
342 for(std::size_t m=0;m<matchedNodes->size();m++) {
343 stk::mesh::EntityId oldNodeId = (*matchedNodes)[m].first;
344 std::size_t newNodeId = (*matchedNodes)[m].second;
346 std::vector<stk::mesh::Entity> elements;
347 std::vector<int> localIds;
351 if((*matchTypes)[m] == 0)
352 offset1 = nodeOffset-offset0;
353 else if((*matchTypes)[m] == 1){
354 offset0 = stkMeshDB_->getMaxEntityId(stkMeshDB_->getNodeRank());
355 offset1 = edgeOffset-offset0;
356 }
else if((*matchTypes)[m] == 2){
357 offset0 = stkMeshDB_->getMaxEntityId(stkMeshDB_->getNodeRank())+stkMeshDB_->getMaxEntityId(stkMeshDB_->getEdgeRank());
358 offset1 = faceOffset-offset0;
363 stkMeshDB_->getOwnedElementsSharingNode(oldNodeId-offset0,elements,localIds,(*matchTypes)[m]);
366 for(std::size_t e=0;e<elements.size();e++){
367 modifySubcellConnectivities(fp,elements[e],(*matchTypes)[m],localIds[e],newNodeId,offset1);
375 template <
typename GO>
378 std::vector<std::size_t> & localCellIds,
379 Kokkos::DynRankView<double,PHX::Device> & points)
const
385 Kokkos::DynRankView<double,PHX::Device> vertices;
389 points = Kokkos::DynRankView<double,PHX::Device>(
"points",localCellIds.size(),numIds,dim);
393 template <
typename GO>
396 return ! sidesetsToAssociate_.empty();
399 template <
typename GO>
402 sidesetsToAssociate_.push_back(sideset_id);
403 sidesetYieldedAssociations_.push_back(
false);
408 stk::mesh::Entity
const e)
410 return static_cast<std::size_t
>(
411 std::distance(elements.begin(), std::find(elements.begin(), elements.end(), e)));
414 template <
typename GO>
417 stk::mesh::BulkData& bulkData = *stkMeshDB_->getBulkData();
418 elmtToAssociatedElmts_.resize(elements_->size());
419 for (std::size_t i = 0; i < sidesetsToAssociate_.size(); ++i) {
420 std::vector<stk::mesh::Entity> sides;
421 stkMeshDB_->getAllSides(sidesetsToAssociate_[i], sides);
422 sidesetYieldedAssociations_[i] = ! sides.empty();
423 for (std::vector<stk::mesh::Entity>::const_iterator si = sides.begin();
424 si != sides.end(); ++si) {
425 stk::mesh::Entity side = *si;
426 const size_t num_elements = bulkData.num_elements(side);
427 stk::mesh::Entity
const* elements = bulkData.begin_elements(side);
428 if (num_elements != 2) {
432 sidesetYieldedAssociations_[i] =
false;
435 const std::size_t ea_id =
getElementIdx(*elements_, elements[0]),
437 elmtToAssociatedElmts_[ea_id].push_back(eb_id);
438 elmtToAssociatedElmts_[eb_id].push_back(ea_id);
443 template <
typename GO>
447 std::vector<std::string> sidesets;
448 for (std::size_t i = 0; i < sidesetYieldedAssociations_.size(); ++i) {
449 int sya, my_sya = sidesetYieldedAssociations_[i] ? 1 : 0;
452 sidesets.push_back(sidesetsToAssociate_[i]);
457 template <
typename GO>
458 const std::vector<typename STKConnManager<GO>::LocalOrdinal>&
461 return elmtToAssociatedElmts_[el];
virtual int getDimension() const =0
panzer::ConnManager< int, GO >::GlobalOrdinal GlobalOrdinal
void associateElementsInSideset(const std::string sideset_id)
virtual Teuchos::RCP< panzer::ConnManagerBase< int > > noConnectivityClone() const
void getIdsAndVertices(const panzer_stk::STK_Interface &mesh, std::string blockId, std::vector< std::size_t > &localIds, ArrayT &vertices)
virtual void buildConnectivity(const panzer::FieldPattern &fp)
void buildLocalElementMapping()
virtual int getDimension() const
STKConnManager(const Teuchos::RCP< const STK_Interface > &stkMeshDB)
virtual bool hasAssociatedNeighbors() const
virtual const std::vector< LocalOrdinal > & getAssociatedNeighbors(const LocalOrdinal &el) const
LocalIdCompare(const RCP< const STK_Interface > &mesh)
void applyInterfaceConditions()
virtual void getDofCoords(const std::string &blockId, const panzer::Intrepid2FieldPattern &coordProvider, std::vector< std::size_t > &localCellIds, Kokkos::DynRankView< double, PHX::Device > &points) const
virtual int numberIds() const
panzer::ConnManager< int, GO >::LocalOrdinal LocalOrdinal
LocalOrdinal addSubcellConnectivities(stk::mesh::Entity element, unsigned subcellRank, LocalOrdinal idCnt, GlobalOrdinal offset)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
void modifySubcellConnectivities(const panzer::FieldPattern &fp, stk::mesh::Entity element, unsigned subcellRank, unsigned subcellId, GlobalOrdinal newId, GlobalOrdinal offset)
std::vector< std::string > checkAssociateElementsInSidesets(const Teuchos::Comm< int > &comm) const
TEUCHOS_DEPRECATED void reduceAll(const Comm< Ordinal > &comm, const EReductionType reductType, const Packet &send, Packet *globalReduct)
void applyPeriodicBCs(const panzer::FieldPattern &fp, GlobalOrdinal nodeOffset, GlobalOrdinal edgeOffset, GlobalOrdinal faceOffset, GlobalOrdinal cellOffset)
virtual std::string getBlockId(LocalOrdinal localElmtId) const
void buildOffsetsAndIdCounts(const panzer::FieldPattern &fp, LocalOrdinal &nodeIdCnt, LocalOrdinal &edgeIdCnt, LocalOrdinal &faceIdCnt, LocalOrdinal &cellIdCnt, GlobalOrdinal &nodeOffset, GlobalOrdinal &edgeOffset, GlobalOrdinal &faceOffset, GlobalOrdinal &cellOffset) const
RCP< const STK_Interface > mesh_
void getInterpolatoryCoordinates(Kokkos::DynRankView< double, PHX::Device > &coords) const
void clearLocalElementMapping()
bool operator()(stk::mesh::Entity a, stk::mesh::Entity b)
virtual const std::vector< int > & getSubcellIndices(int dim, int cellIndex) const =0
#define TEUCHOS_ASSERT(assertion_test)
std::size_t getElementIdx(const std::vector< stk::mesh::Entity > &elements, stk::mesh::Entity const e)