45 #include <stk_mesh/base/GetEntities.hpp>
46 #include <stk_mesh/base/GetBuckets.hpp>
47 #include <stk_mesh/base/FieldBase.hpp>
50 #include "Tpetra_Map.hpp"
51 #include "Tpetra_Import.hpp"
52 #include "Tpetra_Vector.hpp"
54 #include "Teuchos_FancyOStream.hpp"
56 namespace panzer_stk {
58 namespace periodic_helpers {
62 const std::vector<std::pair<std::size_t,std::size_t> > & locallyMatchedIds,
65 using LO = panzer::LocalOrdinal;
66 using GO = panzer::GlobalOrdinal;
68 using Map = Tpetra::Map<LO,GO,NODE>;
69 using Importer = Tpetra::Import<LO,GO,NODE>;
75 int myVal = failure ? 1 : 0;
80 std::vector<GO> requiredInts(locallyRequiredIds.size());
81 for(std::size_t i=0;i<requiredInts.size();i++)
82 requiredInts[i] = locallyRequiredIds[i];
84 std::vector<GO> providedInts(locallyMatchedIds.size());
85 for(std::size_t i=0;i<locallyMatchedIds.size();i++)
86 providedInts[i] = locallyMatchedIds[i].first;
92 Importer importer(providedMap,requiredMap);
95 Tpetra::Vector<GO,LO,GO,NODE> providedVector(providedMap);
97 auto pvHost = providedVector.getLocalViewHost(Tpetra::Access::OverwriteAll);
98 for(std::size_t i=0;i<locallyMatchedIds.size();i++)
99 pvHost(i,0) = locallyMatchedIds[i].second;
102 Tpetra::Vector<GO,LO,GO,NODE> requiredVector(requiredMap);
103 requiredVector.doImport(providedVector,importer,Tpetra::INSERT);
107 =
Teuchos::rcp(
new std::vector<std::pair<std::size_t,std::size_t> >(requiredInts.size()));
109 auto rvHost = requiredVector.getLocalViewHost(Tpetra::Access::ReadOnly);
110 for(std::size_t i=0;i<result->size();i++) {
111 (*result)[i].first = requiredInts[i];
112 (*result)[i].second = rvHost(i,0);
124 const std::string & sideName,
const std::string type_)
131 std::stringstream ss;
132 ss <<
"Can't find part=\"" << sideName <<
"\"" << std::endl;
133 stk::mesh::Part * side = metaData->get_part(sideName,ss.str().c_str());
134 stk::mesh::Selector mySides = *side & (metaData->locally_owned_part() | metaData->globally_shared_part());
136 stk::mesh::EntityRank rank;
137 panzer::GlobalOrdinal offset = 0;
138 if(type_ ==
"coord"){
140 }
else if(type_ ==
"edge"){
143 }
else if(type_ ==
"face"){
147 ss <<
"Can't do BCs of type " << type_ << std::endl;
151 std::vector<stk::mesh::Bucket*>
const& nodeBuckets =
152 bulkData->get_buckets(rank, mySides);
156 std::size_t nodeCount = 0;
157 for(std::size_t b=0;b<nodeBuckets.size();b++)
158 nodeCount += nodeBuckets[b]->size();
161 =
Teuchos::rcp(
new std::vector<std::size_t>(nodeCount));
164 for(std::size_t b=0,index=0;b<nodeBuckets.size();b++) {
165 stk::mesh::Bucket & bucket = *nodeBuckets[b];
167 for(std::size_t n=0;n<bucket.size();n++,index++)
168 (*sideIds)[index] = bulkData->identifier(bucket[n]) + offset;
174 std::pair<Teuchos::RCP<std::vector<std::size_t> >,
177 const std::string & sideName,
const std::string type_)
186 std::stringstream ss;
187 ss <<
"Can't find part=\"" << sideName <<
"\"" << std::endl;
188 stk::mesh::Part * side = metaData->get_part(sideName,ss.str().c_str());
189 stk::mesh::Selector mySides = (*side) & metaData->locally_owned_part();
191 stk::mesh::EntityRank rank;
193 stk::mesh::EntityId offset = 0;
194 if(type_ ==
"coord"){
197 }
else if(type_ ==
"edge"){
201 }
else if(type_ ==
"face"){
206 ss <<
"Can't do BCs of type " << type_ << std::endl;
210 std::vector<stk::mesh::Bucket*>
const& nodeBuckets =
211 bulkData->get_buckets(rank, mySides);
215 std::size_t nodeCount = 0;
216 for(std::size_t b=0;b<nodeBuckets.size();b++)
217 nodeCount += nodeBuckets[b]->size();
220 =
Teuchos::rcp(
new std::vector<std::size_t>(nodeCount));
225 for(std::size_t b=0,index=0;b<nodeBuckets.size();b++) {
226 stk::mesh::Bucket & bucket = *nodeBuckets[b];
227 double const* array = stk::mesh::field_data(*field, bucket);
229 for(std::size_t n=0;n<bucket.size();n++,index++) {
230 (*sideIds)[index] = bulkData->identifier(bucket[n]) + offset;
234 for(std::size_t d=0;d<physicalDim;d++)
235 coord[d] = array[physicalDim*n + d];
239 for(std::size_t d=physicalDim;d<3;d++)
244 return std::make_pair(sideIds,sideCoords);
247 std::pair<Teuchos::RCP<std::vector<std::size_t> >,
250 const std::string & sideName,
const std::string type_)
254 using LO = panzer::LocalOrdinal;
255 using GO = panzer::GlobalOrdinal;
257 using Map = Tpetra::Map<LO,GO,NODE>;
258 using Importer = Tpetra::Import<LO,GO,NODE>;
269 std::pair<Teuchos::RCP<std::vector<std::size_t> >,
273 std::vector<std::size_t> & local_side_ids = *sidePair.first;
274 std::vector<Teuchos::Tuple<double,3> > & local_side_coords = *sidePair.second;
275 std::size_t nodeCount = local_side_ids.size();
279 RCP<Map> idMap_ =
rcp(
new Map(computeInternally,nodeCount,0,comm));
285 auto lidHost = localIdVec_->getLocalViewHost(Tpetra::Access::OverwriteAll);
286 auto lcoordHost = localCoordVec_->getLocalViewHost(Tpetra::Access::OverwriteAll);
287 for(std::size_t n=0;n<local_side_ids.size();n++) {
288 std::size_t nodeId = local_side_ids[n];
291 lidHost(n,0) =
static_cast<GO
>(nodeId);
292 for(
unsigned d=0;d<physicalDim;d++)
293 lcoordHost(n,d) = coords[d];
301 std::size_t dist_nodeCount = idMap_->getGlobalNumElements();
304 RCP<Map> distMap_ =
rcp(
new Map(dist_nodeCount,0,comm,Tpetra::LocallyReplicated));
309 Importer importer_(idMap_,distMap_);
310 distIdVec_->doImport(*localIdVec_,importer_,Tpetra::INSERT);
311 distCoordVec_->doImport(*localCoordVec_,importer_,Tpetra::INSERT);
317 =
Teuchos::rcp(
new std::vector<std::size_t>(dist_nodeCount));
322 const auto didHost = distIdVec_->getLocalViewHost(Tpetra::Access::ReadOnly);
323 const auto dcoordHost = distCoordVec_->getLocalViewHost(Tpetra::Access::ReadOnly);
324 for(std::size_t n=0;n<dist_side_ids->size();++n) {
325 (*dist_side_ids)[n] = didHost(n,0);
328 for(
unsigned d=0;d<physicalDim;++d) {
329 coords[d] = dcoordHost(n,d);
332 for(
unsigned d=physicalDim;d<3;++d)
336 return std::make_pair(dist_side_ids,dist_side_coords);
stk::mesh::Field< double > VectorFieldType
Teuchos::RCP< std::vector< std::size_t > > getLocalSideIds(const STK_Interface &mesh, const std::string &sideName, const std::string type_)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
stk::mesh::EntityId getMaxEntityId(unsigned entityRank) const
get max entity ID of type entityRank
std::pair< Teuchos::RCP< std::vector< std::size_t > >, Teuchos::RCP< std::vector< Teuchos::Tuple< double, 3 > > > > getLocalSideIdsAndCoords(const STK_Interface &mesh, const std::string &sideName, const std::string type_)
PHX::MDField< ScalarT, panzer::Cell, panzer::IP > result
A field that will be used to build up the result of the integral we're performing.
unsigned getDimension() const
get the dimension
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Teuchos::RCP< stk::mesh::BulkData > getBulkData() const
stk::mesh::EntityRank getFaceRank() const
Teuchos::RCP< stk::mesh::MetaData > getMetaData() const
const VectorFieldType & getEdgesField() const
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...
stk::mesh::EntityRank getEdgeRank() const
Teuchos::RCP< std::vector< std::pair< std::size_t, std::size_t > > > getGlobalPairing(const std::vector< std::size_t > &locallyRequiredIds, const std::vector< std::pair< std::size_t, std::size_t > > &locallyMatchedIds, const STK_Interface &mesh, bool failure)
Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
get the comm associated with this mesh
const VectorFieldType & getFacesField() const
#define TEUCHOS_ASSERT(assertion_test)
stk::mesh::EntityRank getNodeRank() const
Tpetra::KokkosCompat::KokkosDeviceWrapperNode< PHX::Device > TpetraNodeType
std::pair< Teuchos::RCP< std::vector< std::size_t > >, Teuchos::RCP< std::vector< Teuchos::Tuple< double, 3 > > > > getSideIdsAndCoords(const STK_Interface &mesh, const std::string &sideName, const std::string type_)
const VectorFieldType & getCoordinatesField() const