43 #ifndef __Panzer_STK_Interface_hpp__ 44 #define __Panzer_STK_Interface_hpp__ 46 #include <Teuchos_RCP.hpp> 47 #include <Teuchos_DefaultMpiComm.hpp> 49 #include <stk_mesh/base/Types.hpp> 50 #include <stk_mesh/base/MetaData.hpp> 51 #include <stk_mesh/base/BulkData.hpp> 52 #include <stk_mesh/base/Field.hpp> 53 #include <stk_mesh/base/FieldBase.hpp> 54 #include <stk_mesh/base/MetaData.hpp> 55 #include <stk_mesh/base/CoordinateSystems.hpp> 57 #include "Kokkos_Core.hpp" 59 #include <Shards_CellTopology.hpp> 60 #include <Shards_CellTopologyData.h> 62 #include <PanzerAdaptersSTK_config.hpp> 63 #include <Kokkos_ViewFactory.hpp> 65 #include <unordered_map> 67 #ifdef PANZER_HAVE_IOSS 68 #include <stk_io/StkMeshIoBroker.hpp> 71 #ifdef PANZER_HAVE_PERCEPT 74 class URP_Heterogeneous_3D;
80 class PeriodicBC_MatcherBase;
87 ElementDescriptor(stk::mesh::EntityId gid,
const std::vector<stk::mesh::EntityId> & nodes);
94 std::vector<stk::mesh::EntityId>
nodes_;
101 Teuchos::RCP<ElementDescriptor>
131 void addElementBlock(
const std::string & name,
const CellTopologyData * ctData);
135 void addSideset(
const std::string & name,
const CellTopologyData * ctData);
143 void addSolutionField(
const std::string & fieldName,
const std::string & blockId);
147 void addCellField(
const std::string & fieldName,
const std::string & blockId);
158 const std::vector<std::string> & coordField,
159 const std::string & dispPrefix);
174 void initialize(stk::ParallelMachine parallelMach,
bool setupIO=
true,
175 const bool buildRefinementSupport =
false);
200 void addNode(stk::mesh::EntityId gid,
const std::vector<double> & coord);
202 void addElement(
const Teuchos::RCP<ElementDescriptor> & ed,stk::mesh::Part * block);
243 std::vector<stk::mesh::EntityId> & subcellIds)
const;
247 void getMyElements(std::vector<stk::mesh::Entity> & elements)
const;
251 void getMyElements(
const std::string & blockID,std::vector<stk::mesh::Entity> & elements)
const;
260 void getNeighborElements(
const std::string & blockID,std::vector<stk::mesh::Entity> & elements)
const;
269 void getMySides(
const std::string & sideName,std::vector<stk::mesh::Entity> & sides)
const;
279 void getMySides(
const std::string & sideName,
const std::string & blockName,std::vector<stk::mesh::Entity> & sides)
const;
288 void getAllSides(
const std::string & sideName,std::vector<stk::mesh::Entity> & sides)
const;
299 void getAllSides(
const std::string & sideName,
const std::string & blockName,std::vector<stk::mesh::Entity> & sides)
const;
310 void getMyNodes(
const std::string & sideName,
const std::string & blockName,std::vector<stk::mesh::Entity> & nodes)
const;
320 stk::mesh::Entity
findConnectivityById(stk::mesh::Entity src, stk::mesh::EntityRank tgt_rank,
unsigned rel_id)
const;
334 const std::string& filename);
352 const std::string& filename);
391 const std::string& key,
415 const std::string& key,
416 const double& value);
439 const std::string& key,
440 const std::vector<int>& value);
463 const std::string& key,
464 const std::vector<double>& value);
470 Teuchos::RCP<const Teuchos::Comm<int> >
getComm()
const;
478 {
if(
bulkData_==Teuchos::null)
return false;
479 return bulkData_->in_modifiable_state(); }
523 std::map<std::string, stk::mesh::Part*>::const_iterator itr =
elementBlocks_.find(name);
535 return (itr !=
sidesets_.end()) ? itr->second :
nullptr;
545 return (itr !=
nodesets_.end()) ? itr->second :
nullptr;
561 void getNodeIdsForElement(stk::mesh::Entity element,std::vector<stk::mesh::EntityId> & nodeIds)
const;
566 std::vector<int> & relIds)
const;
571 std::vector<int> & relIds,
unsigned int matchType)
const;
575 void getElementsSharingNodes(
const std::vector<stk::mesh::EntityId> nodeId,std::vector<stk::mesh::Entity> & elements)
const;
601 {
return bulkData_->parallel_owner_rank(entity); }
605 inline bool isValid(stk::mesh::Entity entity)
const 617 const std::string & blockId)
const;
623 stk::mesh::Field<double> *
getCellField(
const std::string & fieldName,
624 const std::string & blockId)
const;
650 template <
typename ArrayT>
652 const std::vector<std::size_t> & localElementIds,
const ArrayT & solutionValues,
double scaleValue=1.0);
668 template <
typename ArrayT>
670 const std::vector<std::size_t> & localElementIds,ArrayT & solutionValues)
const;
686 template <
typename ArrayT>
687 void setCellFieldData(
const std::string & fieldName,
const std::string & blockId,
688 const std::vector<std::size_t> & localElementIds,
const ArrayT & solutionValues,
double scaleValue=1.0);
698 template <
typename ArrayT>
699 void getElementVertices(
const std::vector<std::size_t> & localIds, ArrayT & vertices)
const;
709 template <
typename ArrayT>
710 void getElementVertices(
const std::vector<stk::mesh::Entity> & elements, ArrayT & vertices)
const;
721 template <
typename ArrayT>
722 void getElementVertices(
const std::vector<std::size_t> & localIds,
const std::string & eBlock, ArrayT & vertices)
const;
733 template <
typename ArrayT>
734 void getElementVertices(
const std::vector<stk::mesh::Entity> & elements,
const std::string & eBlock, ArrayT & vertices)
const;
744 template <
typename ArrayT>
755 template <
typename ArrayT>
767 template <
typename ArrayT>
768 void getElementVerticesNoResize(
const std::vector<std::size_t> & localIds,
const std::string & eBlock, ArrayT & vertices)
const;
779 template <
typename ArrayT>
780 void getElementVerticesNoResize(
const std::vector<stk::mesh::Entity> & elements,
const std::string & eBlock, ArrayT & vertices)
const;
785 stk::mesh::EntityRank
getElementRank()
const {
return stk::topology::ELEMENT_RANK; }
787 stk::mesh::EntityRank
getFaceRank()
const {
return stk::topology::FACE_RANK; }
788 stk::mesh::EntityRank
getEdgeRank()
const {
return stk::topology::EDGE_RANK; }
789 stk::mesh::EntityRank
getNodeRank()
const {
return stk::topology::NODE_RANK; }
801 const std::vector<Teuchos::RCP<const PeriodicBC_MatcherBase> > &
807 std::vector<Teuchos::RCP<const PeriodicBC_MatcherBase> > &
824 void addPeriodicBCs(
const std::vector<Teuchos::RCP<const PeriodicBC_MatcherBase> > & bc_vec)
827 std::pair<Teuchos::RCP<std::vector<std::pair<std::size_t,std::size_t> > >, Teuchos::RCP<std::vector<unsigned int> > >
836 void print(std::ostream & os)
const;
844 Teuchos::RCP<const shards::CellTopology>
getCellTopology(
const std::string & eBlock)
const;
867 void rebalance(
const Teuchos::ParameterList & params);
906 template <
typename ArrayT>
909 template <
typename ArrayT>
911 const std::string & eBlock, ArrayT & vertices)
const;
922 template <
typename ArrayT>
925 template <
typename ArrayT>
933 void refineMesh(
const int numberOfLevels,
const bool deleteParentElements);
961 Teuchos::RCP<Teuchos::MpiComm<int> >
getSafeCommunicator(stk::ParallelMachine parallelMach)
const;
973 bool isMeshCoordField(
const std::string & eBlock,
const std::string & fieldName,
int & axis)
const;
990 template <
typename ArrayT>
991 void setDispFieldData(
const std::string & fieldName,
const std::string & blockId,
int axis,
992 const std::vector<std::size_t> & localElementIds,
const ArrayT & solutionValues);
998 #ifdef PANZER_HAVE_PERCEPT 999 Teuchos::RCP<percept::PerceptMesh> refinedMesh_;
1000 Teuchos::RCP<percept::URP_Heterogeneous_3D> breakPattern_;
1044 #ifdef PANZER_HAVE_IOSS 1046 Teuchos::RCP<stk::io::StkMeshIoBroker> meshData_;
1053 enum class GlobalVariable
1077 const GlobalVariable& flag);
1082 Teuchos::ParameterList globalData_;
1118 template <
typename ArrayT>
1120 const std::vector<std::size_t> & localElementIds,
const ArrayT & solutionValues,
double scaleValue)
1124 int field_axis = -1;
1126 setDispFieldData(fieldName,blockId,field_axis,localElementIds,solutionValues);
1132 for(std::size_t cell=0;cell<localElementIds.size();cell++) {
1133 std::size_t localId = localElementIds[cell];
1134 stk::mesh::Entity element = elements[localId];
1137 const size_t num_nodes =
bulkData_->num_nodes(element);
1138 stk::mesh::Entity
const* nodes =
bulkData_->begin_nodes(element);
1139 for(std::size_t i=0; i<num_nodes; ++i) {
1140 stk::mesh::Entity node = nodes[i];
1142 double * solnData = stk::mesh::field_data(*
field,node);
1144 solnData[0] = scaleValue*solutionValues(cell,i);
1149 template <
typename ArrayT>
1151 const std::vector<std::size_t> & localElementIds,
const ArrayT & dispValues)
1153 TEUCHOS_ASSERT(axis>=0);
1160 for(std::size_t cell=0;cell<localElementIds.size();cell++) {
1161 std::size_t localId = localElementIds[cell];
1162 stk::mesh::Entity element = elements[localId];
1165 const size_t num_nodes =
bulkData_->num_nodes(element);
1166 stk::mesh::Entity
const* nodes =
bulkData_->begin_nodes(element);
1167 for(std::size_t i=0; i<num_nodes; ++i) {
1168 stk::mesh::Entity node = nodes[i];
1170 double * solnData = stk::mesh::field_data(*
field,node);
1171 double * coordData = stk::mesh::field_data(coord_field,node);
1173 solnData[0] = dispValues(cell,i)-coordData[axis];
1178 template <
typename ArrayT>
1180 const std::vector<std::size_t> & localElementIds,ArrayT & solutionValues)
const 1186 localElementIds.size(),
1187 bulkData_->num_nodes(elements[localElementIds[0]]));
1191 for(std::size_t cell=0;cell<localElementIds.size();cell++) {
1192 std::size_t localId = localElementIds[cell];
1193 stk::mesh::Entity element = elements[localId];
1196 const size_t num_nodes =
bulkData_->num_nodes(element);
1197 stk::mesh::Entity
const* nodes =
bulkData_->begin_nodes(element);
1198 for(std::size_t i=0; i<num_nodes; ++i) {
1199 stk::mesh::Entity node = nodes[i];
1201 double * solnData = stk::mesh::field_data(*
field,node);
1203 solutionValues(cell,i) = solnData[0];
1208 template <
typename ArrayT>
1210 const std::vector<std::size_t> & localElementIds,
const ArrayT & solutionValues,
double scaleValue)
1216 for(std::size_t cell=0;cell<localElementIds.size();cell++) {
1217 std::size_t localId = localElementIds[cell];
1218 stk::mesh::Entity element = elements[localId];
1220 double * solnData = stk::mesh::field_data(*
field,element);
1221 TEUCHOS_ASSERT(solnData!=0);
1222 solnData[0] = scaleValue*solutionValues.access(cell,0);
1226 template <
typename ArrayT>
1237 std::vector<stk::mesh::Entity> selected_elements;
1238 for(std::size_t cell=0;cell<localElementIds.size();cell++)
1239 selected_elements.push_back(elements[localElementIds[cell]]);
1244 TEUCHOS_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
1245 "STK_Interface::getElementVertices: Cannot call this method when field coordinates are used " 1246 "without specifying an element block.");
1250 template <
typename ArrayT>
1257 TEUCHOS_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
1258 "STK_Interface::getElementVertices: Cannot call this method when field coordinates are used " 1259 "without specifying an element block.");
1263 template <
typename ArrayT>
1274 template <
typename ArrayT>
1280 std::vector<stk::mesh::Entity> selected_elements;
1281 for(std::size_t cell=0;cell<localElementIds.size();cell++)
1282 selected_elements.push_back(elements[localElementIds[cell]]);
1292 template <
typename ArrayT>
1303 std::vector<stk::mesh::Entity> selected_elements;
1304 for(std::size_t cell=0;cell<localElementIds.size();cell++)
1305 selected_elements.push_back(elements[localElementIds[cell]]);
1310 TEUCHOS_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
1311 "STK_Interface::getElementVerticesNoResize: Cannot call this method when field coordinates are used " 1312 "without specifying an element block.");
1316 template <
typename ArrayT>
1323 TEUCHOS_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
1324 "STK_Interface::getElementVerticesNoResize: Cannot call this method when field coordinates are used " 1325 "without specifying an element block.");
1329 template <
typename ArrayT>
1340 template <
typename ArrayT>
1346 std::vector<stk::mesh::Entity> selected_elements;
1347 for(std::size_t cell=0;cell<localElementIds.size();cell++)
1348 selected_elements.push_back(elements[localElementIds[cell]]);
1358 template <
typename ArrayT>
1362 if(elements.size()==0) {
1372 unsigned masterVertexCount
1373 = stk::mesh::get_cell_topology(
bulkData_->bucket(elements[0]).topology()).getCellTopologyData()->vertex_count;
1380 for(std::size_t cell=0;cell<elements.size();cell++) {
1381 stk::mesh::Entity element = elements[cell];
1382 TEUCHOS_ASSERT(element!=0);
1384 unsigned vertexCount
1385 = stk::mesh::get_cell_topology(
bulkData_->bucket(element).topology()).getCellTopologyData()->vertex_count;
1386 TEUCHOS_TEST_FOR_EXCEPTION(vertexCount!=masterVertexCount,std::runtime_error,
1387 "In call to STK_Interface::getElementVertices all elements " 1388 "must have the same vertex count!");
1391 const size_t num_nodes =
bulkData_->num_nodes(element);
1392 stk::mesh::Entity
const* nodes =
bulkData_->begin_nodes(element);
1393 TEUCHOS_TEST_FOR_EXCEPTION(num_nodes!=masterVertexCount,std::runtime_error,
1394 "In call to STK_Interface::getElementVertices cardinality of " 1395 "element node relations must be the vertex count!");
1396 for(std::size_t node=0; node<num_nodes; ++node) {
1400 for(
unsigned d=0;d<dim;d++)
1401 vertices(cell,node,d) = coord[d];
1406 template <
typename ArrayT>
1410 if(elements.size()==0) {
1419 unsigned masterVertexCount
1420 = stk::mesh::get_cell_topology(
bulkData_->bucket(elements[0]).topology()).getCellTopologyData()->vertex_count;
1424 for(std::size_t cell=0;cell<elements.size();cell++) {
1425 stk::mesh::Entity element = elements[cell];
1426 TEUCHOS_ASSERT(element!=0);
1428 unsigned vertexCount
1429 = stk::mesh::get_cell_topology(
bulkData_->bucket(element).topology()).getCellTopologyData()->vertex_count;
1430 TEUCHOS_TEST_FOR_EXCEPTION(vertexCount!=masterVertexCount,std::runtime_error,
1431 "In call to STK_Interface::getElementVertices all elements " 1432 "must have the same vertex count!");
1435 const size_t num_nodes =
bulkData_->num_nodes(element);
1436 stk::mesh::Entity
const* nodes =
bulkData_->begin_nodes(element);
1437 TEUCHOS_TEST_FOR_EXCEPTION(num_nodes!=masterVertexCount,std::runtime_error,
1438 "In call to STK_Interface::getElementVertices cardinality of " 1439 "element node relations must be the vertex count!");
1440 for(std::size_t node=0; node<num_nodes; ++node) {
1444 for(
unsigned d=0;d<dim;d++)
1445 vertices(cell,node,d) = coord[d];
1450 template <
typename ArrayT>
1456 if(elements.size()==0) {
1462 unsigned masterVertexCount
1463 = stk::mesh::get_cell_topology(
bulkData_->bucket(elements[0]).topology()).getCellTopologyData()->vertex_count;
1468 std::map<std::string,std::vector<std::string> >::const_iterator itr =
meshCoordFields_.find(eBlock);
1471 TEUCHOS_ASSERT(
false);
1474 const std::vector<std::string> & coordField = itr->second;
1475 std::vector<SolutionFieldType*> fields(
getDimension());
1476 for(std::size_t d=0;d<fields.size();d++) {
1480 for(std::size_t cell=0;cell<elements.size();cell++) {
1481 stk::mesh::Entity element = elements[cell];
1484 const size_t num_nodes =
bulkData_->num_nodes(element);
1485 stk::mesh::Entity
const* nodes =
bulkData_->begin_nodes(element);
1486 for(std::size_t i=0; i<num_nodes; ++i) {
1487 stk::mesh::Entity node = nodes[i];
1492 double * solnData = stk::mesh::field_data(*fields[d],node);
1496 vertices(cell,i,d) = solnData[0]+coord[d];
1502 template <
typename ArrayT>
1504 const std::string & eBlock, ArrayT & vertices)
const 1509 if(elements.size()==0) {
1513 std::map<std::string,std::vector<std::string> >::const_iterator itr =
meshCoordFields_.find(eBlock);
1516 TEUCHOS_ASSERT(
false);
1519 const std::vector<std::string> & coordField = itr->second;
1520 std::vector<SolutionFieldType*> fields(
getDimension());
1521 for(std::size_t d=0;d<fields.size();d++) {
1525 for(std::size_t cell=0;cell<elements.size();cell++) {
1526 stk::mesh::Entity element = elements[cell];
1529 const size_t num_nodes =
bulkData_->num_nodes(element);
1530 stk::mesh::Entity
const* nodes =
bulkData_->begin_nodes(element);
1531 for(std::size_t i=0; i<num_nodes; ++i) {
1532 stk::mesh::Entity node = nodes[i];
1537 double * solnData = stk::mesh::field_data(*fields[d],node);
1541 vertices(cell,i,d) = solnData[0]+coord[d];
void initializeFromMetaData()
Teuchos::RCP< ElementDescriptor > buildElementDescriptor(stk::mesh::EntityId elmtId, std::vector< stk::mesh::EntityId > &nodes)
Teuchos::RCP< stk::mesh::MetaData > getMetaData() const
stk::mesh::Part * getElementBlockPart(const std::string &name) const
get the block count
bool isMeshCoordField(const std::string &eBlock, const std::string &fieldName, int &axis) const
ElementBlockException(const std::string &what)
void addNodeset(const std::string &name)
void getElementVertices(const std::vector< std::size_t > &localIds, ArrayT &vertices) const
Teuchos::RCP< stk::mesh::BulkData > bulkData_
void setSolutionFieldData(const std::string &fieldName, const std::string &blockId, const std::vector< std::size_t > &localElementIds, const ArrayT &solutionValues, double scaleValue=1.0)
void getMySides(const std::string &sideName, std::vector< stk::mesh::Entity > &sides) const
std::vector< stk::mesh::Part * > facesPartVec_
stk::mesh::EntityId getGID() const
void getElementVertices_FromCoordsNoResize(const std::vector< stk::mesh::Entity > &elements, ArrayT &vertices) const
std::vector< stk::mesh::EntityId > nodes_
std::map< std::pair< std::string, std::string >, SolutionFieldType * > fieldNameToSolution_
bool getUseLowerCaseForIO() const
std::string containingBlockId(stk::mesh::Entity elmt) const
void addEntityToNodeset(stk::mesh::Entity entity, stk::mesh::Part *nodeset)
stk::mesh::EntityRank getFaceRank() const
void addGlobalToExodus(const std::string &key, const int &value)
Add an int global variable to the information to be written to the Exodus output file.
void setDispFieldData(const std::string &fieldName, const std::string &blockId, int axis, const std::vector< std::size_t > &localElementIds, const ArrayT &solutionValues)
void addPeriodicBC(const Teuchos::RCP< const PeriodicBC_MatcherBase > &bc)
std::size_t elementLocalId(stk::mesh::Entity elmt) const
std::vector< Teuchos::RCP< const PeriodicBC_MatcherBase > > & getPeriodicBCVector()
stk::mesh::Field< double > SolutionFieldType
std::vector< stk::mesh::Part * > edgesPartVec_
const std::vector< stk::mesh::EntityId > & getNodes() const
std::map< std::string, double > blockWeights_
VectorFieldType * coordinatesField_
bool isModifiable() const
void getMyElements(std::vector< stk::mesh::Entity > &elements) const
void addEntityToSideset(stk::mesh::Entity entity, stk::mesh::Part *sideset)
stk::mesh::Part * nodesPart_
std::size_t getNumSidesets() const
get the side set count
stk::mesh::EntityRank getNodeRank() const
void addSolutionField(const std::string &fieldName, const std::string &blockId)
std::map< std::string, Teuchos::RCP< shards::CellTopology > > elementBlockCT_
Teuchos::RCP< stk::mesh::BulkData > getBulkData() const
std::size_t currentLocalId_
bool getUseFieldCoordinates() const
std::map< std::string, stk::mesh::Part * > sidesets_
std::vector< std::size_t > entityCounts_
void getElementsSharingNode(stk::mesh::EntityId nodeId, std::vector< stk::mesh::Entity > &elements) const
get a set of elements sharing a single node
stk::mesh::Part * getOwnedPart() const
Get a pointer to the locally owned part.
void getSubcellIndices(unsigned entityRank, stk::mesh::EntityId elementId, std::vector< stk::mesh::EntityId > &subcellIds) const
stk::mesh::EntityRank getSideRank() const
SolutionFieldType * loadBalField_
stk::mesh::Entity findConnectivityById(stk::mesh::Entity src, stk::mesh::EntityRank tgt_rank, unsigned rel_id) const
bool isValid(stk::mesh::Entity entity) const
static const std::string nodesString
const double * getNodeCoordinates(stk::mesh::EntityId nodeId) const
void instantiateBulkData(stk::ParallelMachine parallelMach)
void getElementBlockNames(std::vector< std::string > &names) const
stk::mesh::EntityRank getEdgeRank() const
stk::mesh::Part * getSideset(const std::string &name) const
unsigned getDimension() const
get the dimension
bool useFieldCoordinates_
const VectorFieldType & getFacesField() const
void addElement(const Teuchos::RCP< ElementDescriptor > &ed, stk::mesh::Part *block)
std::map< std::string, std::vector< std::string > > meshDispFields_
Kokkos::DynRankView< typename InputArray::value_type, PHX::Device > createDynRankView(const InputArray &a, const std::string &name, const DimensionPack... dims)
Wrapper to simplify Panzer use of Sacado ViewFactory.
VectorFieldType * edgesField_
void setUseFieldCoordinates(bool useFieldCoordinates)
stk::mesh::EntityId elementGlobalId(std::size_t lid) const
void getElementVertices_FromField(const std::vector< stk::mesh::Entity > &elements, const std::string &eBlock, ArrayT &vertices) const
void getNodeIdsForElement(stk::mesh::Entity element, std::vector< stk::mesh::EntityId > &nodeIds) const
get a list of node ids for nodes connected to an element
void initialize(stk::ParallelMachine parallelMach, bool setupIO=true, const bool buildRefinementSupport=false)
static const std::string coordsString
stk::mesh::Field< double, stk::mesh::Cartesian > VectorFieldType
stk::mesh::Field< double > * getCellField(const std::string &fieldName, const std::string &blockId) const
std::size_t getEntityCounts(unsigned entityRank) const
get the global counts for the entity of specified rank
void getElementVertices_FromCoords(const std::vector< stk::mesh::Entity > &elements, ArrayT &vertices) const
void writeToExodus(const std::string &filename)
Write this mesh and associated fields to the given output file.
void getElementVertices_FromFieldNoResize(const std::vector< stk::mesh::Entity > &elements, const std::string &eBlock, ArrayT &vertices) const
std::map< std::string, stk::mesh::Part * > elementBlocks_
bool isInitialized() const
Has initialize been called on this mesh object?
std::vector< stk::mesh::Part * > nodesPartVec_
static const std::string edgesString
std::size_t getNumElementBlocks() const
get the block count
SidesetException(const std::string &what)
void printMetaData(std::ostream &os) const
void initializeFieldsInSTK(const std::map< std::pair< std::string, std::string >, SolutionFieldType *> &nameToField, bool setupIO)
double getInitialStateTime() const
const STK_Interface * mesh_
std::vector< Teuchos::RCP< const PeriodicBC_MatcherBase > > periodicBCs_
Teuchos::RCP< const std::vector< stk::mesh::Entity > > getElementsOrderedByLID() const
void addNode(stk::mesh::EntityId gid, const std::vector< double > &coord)
void setBlockWeight(const std::string &blockId, double weight)
const std::vector< Teuchos::RCP< const PeriodicBC_MatcherBase > > & getPeriodicBCVector() const
void buildSubcells()
force the mesh to build subcells: edges and faces
Teuchos::RCP< Teuchos::MpiComm< int > > mpiComm_
void getElementVerticesNoResize(const std::vector< std::size_t > &localIds, ArrayT &vertices) const
LocalIdCompare(const STK_Interface *mesh)
void setInitialStateTime(double value)
std::size_t getNumNodesets() const
get the side set count
Teuchos::RCP< stk::mesh::MetaData > metaData_
virtual ~ElementDescriptor()
Teuchos::RCP< std::vector< stk::mesh::Entity > > orderedElementVector_
Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
get the comm associated with this mesh
void addMeshCoordFields(const std::string &blockId, const std::vector< std::string > &coordField, const std::string &dispPrefix)
void applyElementLoadBalanceWeights()
void setUseLowerCaseForIO(bool useLowerCase)
bool validBlockId(const std::string &blockId) const
void addSideset(const std::string &name, const CellTopologyData *ctData)
void getSolutionFieldData(const std::string &fieldName, const std::string &blockId, const std::vector< std::size_t > &localElementIds, ArrayT &solutionValues) const
void setCellFieldData(const std::string &fieldName, const std::string &blockId, const std::vector< std::size_t > &localElementIds, const ArrayT &solutionValues, double scaleValue=1.0)
ProcIdFieldType * getProcessorIdField()
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...
void rebalance(const Teuchos::ParameterList ¶ms)
Teuchos::RCP< const shards::CellTopology > getCellTopology(const std::string &eBlock) const
void getAllSides(const std::string &sideName, std::vector< stk::mesh::Entity > &sides) const
void getSidesetNames(std::vector< std::string > &name) const
std::map< std::string, std::vector< std::string > > meshCoordFields_
static const std::string facesString
Teuchos::RCP< Teuchos::MpiComm< int > > getSafeCommunicator(stk::ParallelMachine parallelMach) const
std::pair< Teuchos::RCP< std::vector< std::pair< std::size_t, std::size_t > > >, Teuchos::RCP< std::vector< unsigned int > > > getPeriodicNodePairing() const
stk::mesh::EntityId elementGlobalId(stk::mesh::Entity elmt) const
void refineMesh(const int numberOfLevels, const bool deleteParentElements)
void getNodesetNames(std::vector< std::string > &name) const
void buildLocalElementIDs()
void setupExodusFile(const std::string &filename)
Set up an output Exodus file for writing results.
std::unordered_map< stk::mesh::EntityId, std::size_t > localIDHash_
std::vector< stk::mesh::EntityId > maxEntityId_
unsigned entityOwnerRank(stk::mesh::Entity entity) const
const VectorFieldType & getCoordinatesField() const
void print(std::ostream &os) const
stk::mesh::EntityRank getElementRank() const
stk::mesh::Field< ProcIdData > ProcIdFieldType
void addPeriodicBCs(const std::vector< Teuchos::RCP< const PeriodicBC_MatcherBase > > &bc_vec)
stk::mesh::EntityId getMaxEntityId(unsigned entityRank) const
get max entity ID of type entityRank
void addCellField(const std::string &fieldName, const std::string &blockId)
stk::mesh::Part * edgesPart_
bool operator()(stk::mesh::Entity a, stk::mesh::Entity b)
void getMyNodes(const std::string &sideName, const std::string &blockName, std::vector< stk::mesh::Entity > &nodes) const
void getElementsSharingNodes(const std::vector< stk::mesh::EntityId > nodeId, std::vector< stk::mesh::Entity > &elements) const
get a set of elements sharing multiple nodes
stk::mesh::Part * facesPart_
ProcIdFieldType * processorIdField_
void addElementBlock(const std::string &name, const CellTopologyData *ctData)
void getOwnedElementsSharingNode(stk::mesh::Entity node, std::vector< stk::mesh::Entity > &elements, std::vector< int > &relIds) const
double getCurrentStateTime() const
stk::mesh::Part * getNodeset(const std::string &name) const
stk::mesh::Field< double > * getSolutionField(const std::string &fieldName, const std::string &blockId) const
std::map< std::pair< std::string, std::string >, SolutionFieldType * > fieldNameToCellField_
const VectorFieldType & getEdgesField() const
std::map< std::string, stk::mesh::Part * > nodesets_
VectorFieldType * facesField_
void getNeighborElements(std::vector< stk::mesh::Entity > &elements) const