50 #include "Kokkos_DynRankView.hpp" 57 #include "Teuchos_FancyOStream.hpp" 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());
196 edgeOffset = nodeOffset+(maxNodeId+1)*nodeIdCnt;
197 faceOffset = edgeOffset+(maxEdgeId+1)*edgeIdCnt;
198 cellOffset = faceOffset+(maxFaceId+1)*faceIdCnt;
202 && edgeOffset <= faceOffset
203 && faceOffset <= cellOffset);
206 template <
typename GO>
215 stk::mesh::BulkData& bulkData = *stkMeshDB_->getBulkData();
216 const stk::mesh::EntityRank rank =
static_cast<stk::mesh::EntityRank
>(subcellRank);
217 const size_t num_rels = bulkData.num_connectivity(element, rank);
218 stk::mesh::Entity
const* relations = bulkData.begin(element, rank);
219 for(std::size_t sc=0; sc<num_rels; ++sc) {
220 stk::mesh::Entity subcell = relations[sc];
224 connectivity_.push_back(
offset+idCnt*(bulkData.identifier(subcell)-1)+i);
231 template <
typename GO>
236 LocalOrdinal elmtLID = stkMeshDB_->elementLocalId(element);
238 const std::vector<int> & subCellIndices = fp.
getSubcellIndices(subcellRank,subcellId);
241 for(std::size_t i=0;i<subCellIndices.size();i++) {
242 conn[subCellIndices[i]] =
offset+subCellIndices.size()*(newId-1)+i;
246 template <
typename GO>
249 stk::mesh::BulkData& bulkData = *stkMeshDB_->getBulkData();
253 buildLocalElementMapping();
259 LocalOrdinal nodeIdCnt=0, edgeIdCnt=0, faceIdCnt=0, cellIdCnt=0;
260 GlobalOrdinal nodeOffset=0, edgeOffset=0, faceOffset=0, cellOffset=0;
261 buildOffsetsAndIdCounts(fp, nodeIdCnt, edgeIdCnt, faceIdCnt, cellIdCnt,
262 nodeOffset, edgeOffset, faceOffset, cellOffset);
270 for(std::size_t elmtLid=0;elmtLid!=elements_->size();++elmtLid) {
272 stk::mesh::Entity element = (*elements_)[elmtLid];
275 elmtLidToConn_[elmtLid] = connectivity_.size();
278 numIds += addSubcellConnectivities(element,stkMeshDB_->getNodeRank(),nodeIdCnt,nodeOffset);
279 numIds += addSubcellConnectivities(element,stkMeshDB_->getEdgeRank(),edgeIdCnt,edgeOffset);
280 numIds += addSubcellConnectivities(element,stkMeshDB_->getFaceRank(),faceIdCnt,faceOffset);
286 connectivity_.push_back(cellOffset+cellIdCnt*(bulkData.identifier(element)-1));
291 connSize_[elmtLid] = numIds;
294 applyPeriodicBCs( fp, nodeOffset, edgeOffset, faceOffset, cellOffset);
299 if (hasAssociatedNeighbors())
300 applyInterfaceConditions();
303 template <
typename GO>
307 stk::mesh::Entity element = (*elements_)[localElmtId];
309 return stkMeshDB_->containingBlockId(element);
312 template <
typename GO>
320 = stkMeshDB_->getPeriodicNodePairing();
323 = matchedValues.first;
325 = matchedValues.second;
328 if(matchedNodes==Teuchos::null)
return;
330 for(std::size_t m=0;m<matchedNodes->size();m++) {
331 stk::mesh::EntityId oldNodeId = (*matchedNodes)[m].first;
332 std::size_t newNodeId = (*matchedNodes)[m].second;
334 std::vector<stk::mesh::Entity> elements;
335 std::vector<int> localIds;
339 if((*matchTypes)[m] == 0)
340 offset1 = nodeOffset-offset0;
341 else if((*matchTypes)[m] == 1){
342 offset0 = stkMeshDB_->getMaxEntityId(stkMeshDB_->getNodeRank());
343 offset1 = edgeOffset-offset0;
344 }
else if((*matchTypes)[m] == 2){
345 offset0 = stkMeshDB_->getMaxEntityId(stkMeshDB_->getNodeRank())+stkMeshDB_->getMaxEntityId(stkMeshDB_->getEdgeRank());
346 offset1 = faceOffset-offset0;
351 stkMeshDB_->getOwnedElementsSharingNode(oldNodeId-offset0,elements,localIds,(*matchTypes)[m]);
354 for(std::size_t e=0;e<elements.size();e++){
355 modifySubcellConnectivities(fp,elements[e],(*matchTypes)[m],localIds[e],newNodeId,offset1);
363 template <
typename GO>
366 std::vector<std::size_t> & localCellIds,
367 Kokkos::DynRankView<double,PHX::Device> & points)
const 373 Kokkos::DynRankView<double,PHX::Device> vertices;
377 points = Kokkos::DynRankView<double,PHX::Device>(
"points",localCellIds.size(),numIds,dim);
381 template <
typename GO>
384 return ! sidesetsToAssociate_.empty();
387 template <
typename GO>
390 sidesetsToAssociate_.push_back(sideset_id);
391 sidesetYieldedAssociations_.push_back(
false);
396 stk::mesh::Entity
const e)
398 return static_cast<std::size_t
>(
399 std::distance(elements.begin(), std::find(elements.begin(), elements.end(), e)));
402 template <
typename GO>
405 stk::mesh::BulkData& bulkData = *stkMeshDB_->getBulkData();
406 elmtToAssociatedElmts_.resize(elements_->size());
407 for (std::size_t i = 0; i < sidesetsToAssociate_.size(); ++i) {
408 std::vector<stk::mesh::Entity> sides;
409 stkMeshDB_->getAllSides(sidesetsToAssociate_[i], sides);
410 sidesetYieldedAssociations_[i] = ! sides.empty();
411 for (std::vector<stk::mesh::Entity>::const_iterator si = sides.begin();
412 si != sides.end(); ++si) {
413 stk::mesh::Entity side = *si;
414 const size_t num_elements = bulkData.num_elements(side);
415 stk::mesh::Entity
const* elements = bulkData.begin_elements(side);
416 if (num_elements != 2) {
420 sidesetYieldedAssociations_[i] =
false;
423 const std::size_t ea_id =
getElementIdx(*elements_, elements[0]),
425 elmtToAssociatedElmts_[ea_id].push_back(eb_id);
426 elmtToAssociatedElmts_[eb_id].push_back(ea_id);
431 template <
typename GO>
435 std::vector<std::string> sidesets;
436 for (std::size_t i = 0; i < sidesetYieldedAssociations_.size(); ++i) {
437 int sya, my_sya = sidesetYieldedAssociations_[i] ? 1 : 0;
440 sidesets.push_back(sidesetsToAssociate_[i]);
445 template <
typename GO>
446 const std::vector<typename STKConnManager<GO>::LocalOrdinal>&
449 return elmtToAssociatedElmts_[el];
virtual int getDimension() const =0
virtual const std::vector< LocalOrdinal > & getAssociatedNeighbors(const LocalOrdinal &el) const
panzer::ConnManager< int, GO >::GlobalOrdinal GlobalOrdinal
void associateElementsInSideset(const std::string sideset_id)
void getIdsAndVertices(const panzer_stk::STK_Interface &mesh, std::string blockId, std::vector< std::size_t > &localIds, ArrayT &vertices)
virtual std::string getBlockId(LocalOrdinal localElmtId) const
virtual void buildConnectivity(const panzer::FieldPattern &fp)
void buildLocalElementMapping()
virtual Teuchos::RCP< panzer::ConnManagerBase< int > > noConnectivityClone() const
virtual int getDimension() const
virtual void getDofCoords(const std::string &blockId, const panzer::Intrepid2FieldPattern &coordProvider, std::vector< std::size_t > &localCellIds, Kokkos::DynRankView< double, PHX::Device > &points) const
void buildOffsetsAndIdCounts(const panzer::FieldPattern &fp, LocalOrdinal &nodeIdCnt, LocalOrdinal &edgeIdCnt, LocalOrdinal &faceIdCnt, LocalOrdinal &cellIdCnt, GlobalOrdinal &nodeOffset, GlobalOrdinal &edgeOffset, GlobalOrdinal &faceOffset, GlobalOrdinal &cellOffset) const
LocalIdCompare(const RCP< const STK_Interface > &mesh)
void applyInterfaceConditions()
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)
STKConnManager(const Teuchos::RCP< STK_Interface > &stkMeshDB)
void modifySubcellConnectivities(const panzer::FieldPattern &fp, stk::mesh::Entity element, unsigned subcellRank, unsigned subcellId, GlobalOrdinal newId, GlobalOrdinal offset)
virtual bool hasAssociatedNeighbors() const
void getInterpolatoryCoordinates(Kokkos::DynRankView< double, PHX::Device > &coords) 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)
Teuchos::RCP< const Teuchos::Comm< int > > comm
RCP< const STK_Interface > mesh_
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)
std::vector< std::string > checkAssociateElementsInSidesets(const Teuchos::Comm< int > &comm) const