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>
137 void addElementBlock(
const std::string & name,
const CellTopologyData * ctData);
141 void addEdgeBlock(
const std::string & name,
const CellTopologyData * ctData);
145 void addFaceBlock(
const std::string & name,
const CellTopologyData * ctData);
149 void addSideset(
const std::string & name,
const CellTopologyData * ctData);
157 void addSolutionField(
const std::string & fieldName,
const std::string & blockId);
161 void addCellField(
const std::string & fieldName,
const std::string & blockId);
165 void addEdgeField(
const std::string & fieldName,
const std::string & blockId);
169 void addFaceField(
const std::string & fieldName,
const std::string & blockId);
180 const std::vector<std::string> & coordField,
181 const std::string & dispPrefix);
202 void initialize(stk::ParallelMachine parallelMach,
bool setupIO=
true,
203 const bool buildRefinementSupport =
false);
228 void addNode(stk::mesh::EntityId gid,
const std::vector<double> & coord);
230 void addElement(
const Teuchos::RCP<ElementDescriptor> & ed,stk::mesh::Part * block);
285 std::vector<stk::mesh::EntityId> & subcellIds)
const;
289 void getMyElements(std::vector<stk::mesh::Entity> & elements)
const;
293 void getMyElements(
const std::string & blockID,std::vector<stk::mesh::Entity> & elements)
const;
302 void getNeighborElements(
const std::string & blockID,std::vector<stk::mesh::Entity> & elements)
const;
306 void getMyEdges(std::vector<stk::mesh::Entity> & edges)
const;
315 void getMyEdges(
const std::string & edgeBlockName,std::vector<stk::mesh::Entity> & edges)
const;
325 void getMyEdges(
const std::string & edgeBlockName,
const std::string & blockName,std::vector<stk::mesh::Entity> & edges)
const;
334 void getAllEdges(
const std::string & edgeBlockName,std::vector<stk::mesh::Entity> & edges)
const;
344 void getAllEdges(
const std::string & edgeBlockName,
const std::string & blockName,std::vector<stk::mesh::Entity> & edges)
const;
348 void getMyFaces(std::vector<stk::mesh::Entity> & faces)
const;
357 void getMyFaces(
const std::string & faceBlockName,std::vector<stk::mesh::Entity> & faces)
const;
367 void getMyFaces(
const std::string & faceBlockName,
const std::string & blockName,std::vector<stk::mesh::Entity> & faces)
const;
376 void getAllFaces(
const std::string & faceBlockName,std::vector<stk::mesh::Entity> & faces)
const;
386 void getAllFaces(
const std::string & faceBlockName,
const std::string & blockName,std::vector<stk::mesh::Entity> & faces)
const;
395 void getMySides(
const std::string & sideName,std::vector<stk::mesh::Entity> & sides)
const;
405 void getMySides(
const std::string & sideName,
const std::string & blockName,std::vector<stk::mesh::Entity> & sides)
const;
414 void getAllSides(
const std::string & sideName,std::vector<stk::mesh::Entity> & sides)
const;
424 void getAllSides(
const std::string & sideName,
const std::string & blockName,std::vector<stk::mesh::Entity> & sides)
const;
435 void getMyNodes(
const std::string & sideName,
const std::string & blockName,std::vector<stk::mesh::Entity> & nodes)
const;
445 stk::mesh::Entity
findConnectivityById(stk::mesh::Entity src, stk::mesh::EntityRank tgt_rank,
unsigned rel_id)
const;
460 const bool append =
false);
481 const bool append =
false,
482 const bool append_after_restart_time =
false,
483 const double restart_time = 0.0);
522 const std::string& key,
546 const std::string& key,
547 const double& value);
570 const std::string& key,
571 const std::vector<int>& value);
594 const std::string& key,
595 const std::vector<double>& value);
601 Teuchos::RCP<const Teuchos::Comm<int> >
getComm()
const;
606 #ifdef PANZER_HAVE_PERCEPT 607 Teuchos::RCP<percept::PerceptMesh> getRefinedMesh()
const 609 { TEUCHOS_ASSERT(Teuchos::nonnull(refinedMesh_));
return refinedMesh_; }
615 {
if(
bulkData_==Teuchos::null)
return false;
616 return bulkData_->in_modifiable_state(); }
678 std::map<std::string, stk::mesh::Part*>::const_iterator itr =
elementBlocks_.find(name);
686 std::map<std::string, stk::mesh::Part*>::const_iterator itr =
edgeBlocks_.find(name);
694 std::map<std::string, stk::mesh::Part*>::const_iterator itr =
faceBlocks_.find(name);
706 return (itr !=
sidesets_.end()) ? itr->second :
nullptr;
716 return (itr !=
nodesets_.end()) ? itr->second :
nullptr;
732 void getNodeIdsForElement(stk::mesh::Entity element,std::vector<stk::mesh::EntityId> & nodeIds)
const;
737 std::vector<int> & relIds)
const;
742 std::vector<int> & relIds,
unsigned int matchType)
const;
746 void getElementsSharingNodes(
const std::vector<stk::mesh::EntityId> nodeId,std::vector<stk::mesh::Entity> & elements)
const;
779 std::size_t
edgeLocalId(stk::mesh::Entity elmt)
const;
783 std::size_t
edgeLocalId(stk::mesh::EntityId gid)
const;
805 std::size_t
faceLocalId(stk::mesh::Entity elmt)
const;
809 std::size_t
faceLocalId(stk::mesh::EntityId gid)
const;
824 {
return bulkData_->parallel_owner_rank(entity); }
828 inline bool isValid(stk::mesh::Entity entity)
const 840 const std::string & blockId)
const;
846 stk::mesh::Field<double> *
getCellField(
const std::string & fieldName,
847 const std::string & blockId)
const;
853 stk::mesh::Field<double> *
getEdgeField(
const std::string & fieldName,
854 const std::string & blockId)
const;
860 stk::mesh::Field<double> *
getFaceField(
const std::string & fieldName,
861 const std::string & blockId)
const;
887 template <
typename ArrayT>
889 const std::vector<std::size_t> & localElementIds,
const ArrayT & solutionValues,
double scaleValue=1.0);
905 template <
typename ArrayT>
907 const std::vector<std::size_t> & localElementIds,ArrayT & solutionValues)
const;
923 template <
typename ArrayT>
924 void setCellFieldData(
const std::string & fieldName,
const std::string & blockId,
925 const std::vector<std::size_t> & localElementIds,
const ArrayT & solutionValues,
double scaleValue=1.0);
951 template <
typename ArrayT>
952 void setEdgeFieldData(
const std::string & fieldName,
const std::string & blockId,
953 const std::vector<std::size_t> & localEdgeIds,
const ArrayT & edgeValues,
double scaleValue=1.0);
969 template <
typename ArrayT>
970 void setFaceFieldData(
const std::string & fieldName,
const std::string & blockId,
971 const std::vector<std::size_t> & localFaceIds,
const ArrayT & faceValues,
double scaleValue=1.0);
981 template <
typename ArrayT>
982 void getElementVertices(
const std::vector<std::size_t> & localIds, ArrayT & vertices)
const;
992 template <
typename ArrayT>
993 void getElementVertices(
const std::vector<stk::mesh::Entity> & elements, ArrayT & vertices)
const;
1004 template <
typename ArrayT>
1005 void getElementVertices(
const std::vector<std::size_t> & localIds,
const std::string & eBlock, ArrayT & vertices)
const;
1016 template <
typename ArrayT>
1017 void getElementVertices(
const std::vector<stk::mesh::Entity> & elements,
const std::string & eBlock, ArrayT & vertices)
const;
1027 template <
typename ArrayT>
1038 template <
typename ArrayT>
1050 template <
typename ArrayT>
1051 void getElementVerticesNoResize(
const std::vector<std::size_t> & localIds,
const std::string & eBlock, ArrayT & vertices)
const;
1062 template <
typename ArrayT>
1063 void getElementVerticesNoResize(
const std::vector<stk::mesh::Entity> & elements,
const std::string & eBlock, ArrayT & vertices)
const;
1070 stk::mesh::EntityRank
getFaceRank()
const {
return stk::topology::FACE_RANK; }
1071 stk::mesh::EntityRank
getEdgeRank()
const {
return stk::topology::EDGE_RANK; }
1072 stk::mesh::EntityRank
getNodeRank()
const {
return stk::topology::NODE_RANK; }
1092 const std::vector<Teuchos::RCP<const PeriodicBC_MatcherBase> > &
1098 std::vector<Teuchos::RCP<const PeriodicBC_MatcherBase> > &
1115 void addPeriodicBCs(
const std::vector<Teuchos::RCP<const PeriodicBC_MatcherBase> > & bc_vec)
1118 std::pair<Teuchos::RCP<std::vector<std::pair<std::size_t,std::size_t> > >, Teuchos::RCP<std::vector<unsigned int> > >
1127 void print(std::ostream & os)
const;
1135 Teuchos::RCP<const shards::CellTopology>
getCellTopology(
const std::string & eBlock)
const;
1158 void rebalance(
const Teuchos::ParameterList & params);
1197 template <
typename ArrayT>
1198 void getElementVertices_FromField(
const std::vector<stk::mesh::Entity> & elements,
const std::string & eBlock, ArrayT & vertices)
const;
1200 template <
typename ArrayT>
1202 const std::string & eBlock, ArrayT & vertices)
const;
1213 template <
typename ArrayT>
1216 template <
typename ArrayT>
1224 void refineMesh(
const int numberOfLevels,
const bool deleteParentElements);
1254 Teuchos::RCP<Teuchos::MpiComm<int> >
getSafeCommunicator(stk::ParallelMachine parallelMach)
const;
1266 bool isMeshCoordField(
const std::string & eBlock,
const std::string & fieldName,
int & axis)
const;
1283 template <
typename ArrayT>
1284 void setDispFieldData(
const std::string & fieldName,
const std::string & blockId,
int axis,
1285 const std::vector<std::size_t> & localElementIds,
const ArrayT & solutionValues);
1291 #ifdef PANZER_HAVE_PERCEPT 1292 Teuchos::RCP<percept::PerceptMesh> refinedMesh_;
1293 Teuchos::RCP<percept::URP_Heterogeneous_3D> breakPattern_;
1347 #ifdef PANZER_HAVE_IOSS 1349 Teuchos::RCP<stk::io::StkMeshIoBroker> meshData_;
1356 enum class GlobalVariable
1380 const GlobalVariable& flag);
1385 Teuchos::ParameterList globalData_;
1429 template <
typename ArrayT>
1431 const std::vector<std::size_t> & localElementIds,
const ArrayT & solutionValues,
double scaleValue)
1434 auto solutionValues_h = Kokkos::create_mirror_view(solutionValues);
1435 Kokkos::deep_copy(solutionValues_h, solutionValues);
1437 int field_axis = -1;
1439 setDispFieldData(fieldName,blockId,field_axis,localElementIds,solutionValues_h);
1445 for(std::size_t cell=0;cell<localElementIds.size();cell++) {
1446 std::size_t localId = localElementIds[cell];
1447 stk::mesh::Entity element = elements[localId];
1450 const size_t num_nodes =
bulkData_->num_nodes(element);
1451 stk::mesh::Entity
const* nodes =
bulkData_->begin_nodes(element);
1452 for(std::size_t i=0; i<num_nodes; ++i) {
1453 stk::mesh::Entity node = nodes[i];
1455 double * solnData = stk::mesh::field_data(*
field,node);
1457 solnData[0] = scaleValue*solutionValues_h(cell,i);
1462 template <
typename ArrayT>
1464 const std::vector<std::size_t> & localElementIds,
const ArrayT & dispValues)
1466 TEUCHOS_ASSERT(axis>=0);
1473 for(std::size_t cell=0;cell<localElementIds.size();cell++) {
1474 std::size_t localId = localElementIds[cell];
1475 stk::mesh::Entity element = elements[localId];
1478 const size_t num_nodes =
bulkData_->num_nodes(element);
1479 stk::mesh::Entity
const* nodes =
bulkData_->begin_nodes(element);
1480 for(std::size_t i=0; i<num_nodes; ++i) {
1481 stk::mesh::Entity node = nodes[i];
1483 double * solnData = stk::mesh::field_data(*
field,node);
1484 double * coordData = stk::mesh::field_data(coord_field,node);
1486 solnData[0] = dispValues(cell,i)-coordData[axis];
1491 template <
typename ArrayT>
1493 const std::vector<std::size_t> & localElementIds,ArrayT & solutionValues)
const 1499 localElementIds.size(),
1500 bulkData_->num_nodes(elements[localElementIds[0]]));
1504 for(std::size_t cell=0;cell<localElementIds.size();cell++) {
1505 std::size_t localId = localElementIds[cell];
1506 stk::mesh::Entity element = elements[localId];
1509 const size_t num_nodes =
bulkData_->num_nodes(element);
1510 stk::mesh::Entity
const* nodes =
bulkData_->begin_nodes(element);
1511 for(std::size_t i=0; i<num_nodes; ++i) {
1512 stk::mesh::Entity node = nodes[i];
1514 double * solnData = stk::mesh::field_data(*
field,node);
1516 solutionValues(cell,i) = solnData[0];
1521 template <
typename ArrayT>
1523 const std::vector<std::size_t> & localElementIds,
const ArrayT & solutionValues,
double scaleValue)
1529 auto solutionValues_h = Kokkos::create_mirror_view(solutionValues);
1530 Kokkos::deep_copy(solutionValues_h, solutionValues);
1532 for(std::size_t cell=0;cell<localElementIds.size();cell++) {
1533 std::size_t localId = localElementIds[cell];
1534 stk::mesh::Entity element = elements[localId];
1536 double * solnData = stk::mesh::field_data(*
field,element);
1537 TEUCHOS_ASSERT(solnData!=0);
1538 solnData[0] = scaleValue*solutionValues_h.access(cell,0);
1542 template <
typename ArrayT>
1544 const std::vector<std::size_t> & localEdgeIds,
const ArrayT & edgeValues,
double scaleValue)
1550 auto edgeValues_h = Kokkos::create_mirror_view(edgeValues);
1551 Kokkos::deep_copy(edgeValues_h, edgeValues);
1553 for(std::size_t idx=0;idx<localEdgeIds.size();idx++) {
1554 std::size_t localId = localEdgeIds[idx];
1555 stk::mesh::Entity edge = edges[localId];
1557 double * solnData = stk::mesh::field_data(*
field,edge);
1558 TEUCHOS_ASSERT(solnData!=0);
1559 solnData[0] = scaleValue*edgeValues_h.access(idx,0);
1563 template <
typename ArrayT>
1565 const std::vector<std::size_t> & localFaceIds,
const ArrayT & faceValues,
double scaleValue)
1571 auto faceValues_h = Kokkos::create_mirror_view(faceValues);
1572 Kokkos::deep_copy(faceValues_h, faceValues);
1574 for(std::size_t idx=0;idx<localFaceIds.size();idx++) {
1575 std::size_t localId = localFaceIds[idx];
1576 stk::mesh::Entity face = faces[localId];
1578 double * solnData = stk::mesh::field_data(*
field,face);
1579 TEUCHOS_ASSERT(solnData!=0);
1580 solnData[0] = scaleValue*faceValues_h.access(idx,0);
1584 template <
typename ArrayT>
1595 std::vector<stk::mesh::Entity> selected_elements;
1596 for(std::size_t cell=0;cell<localElementIds.size();cell++)
1597 selected_elements.push_back(elements[localElementIds[cell]]);
1602 TEUCHOS_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
1603 "STK_Interface::getElementVertices: Cannot call this method when field coordinates are used " 1604 "without specifying an element block.");
1608 template <
typename ArrayT>
1615 TEUCHOS_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
1616 "STK_Interface::getElementVertices: Cannot call this method when field coordinates are used " 1617 "without specifying an element block.");
1621 template <
typename ArrayT>
1632 template <
typename ArrayT>
1638 std::vector<stk::mesh::Entity> selected_elements;
1639 for(std::size_t cell=0;cell<localElementIds.size();cell++)
1640 selected_elements.push_back(elements[localElementIds[cell]]);
1650 template <
typename ArrayT>
1661 std::vector<stk::mesh::Entity> selected_elements;
1662 for(std::size_t cell=0;cell<localElementIds.size();cell++)
1663 selected_elements.push_back(elements[localElementIds[cell]]);
1668 TEUCHOS_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
1669 "STK_Interface::getElementVerticesNoResize: Cannot call this method when field coordinates are used " 1670 "without specifying an element block.");
1674 template <
typename ArrayT>
1681 TEUCHOS_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
1682 "STK_Interface::getElementVerticesNoResize: Cannot call this method when field coordinates are used " 1683 "without specifying an element block.");
1687 template <
typename ArrayT>
1698 template <
typename ArrayT>
1704 std::vector<stk::mesh::Entity> selected_elements;
1705 for(std::size_t cell=0;cell<localElementIds.size();cell++)
1706 selected_elements.push_back(elements[localElementIds[cell]]);
1716 template <
typename ArrayT>
1720 if(elements.size() == 0) {
1730 const auto masterVertexCount
1731 = stk::mesh::get_cell_topology(
bulkData_->bucket(elements[0]).topology()).getCellTopologyData()->vertex_count;
1735 auto vertices_h = Kokkos::create_mirror_view(vertices);
1736 Kokkos::deep_copy(vertices_h, vertices);
1740 for(std::size_t cell = 0; cell < elements.size(); cell++) {
1741 const auto element = elements[cell];
1742 TEUCHOS_ASSERT(element != 0);
1744 const auto vertexCount
1745 = stk::mesh::get_cell_topology(
bulkData_->bucket(element).topology()).getCellTopologyData()->vertex_count;
1746 TEUCHOS_TEST_FOR_EXCEPTION(vertexCount != masterVertexCount, std::runtime_error,
1747 "In call to STK_Interface::getElementVertices all elements " 1748 "must have the same vertex count!");
1751 const size_t num_nodes =
bulkData_->num_nodes(element);
1752 auto const* nodes =
bulkData_->begin_nodes(element);
1753 TEUCHOS_TEST_FOR_EXCEPTION(num_nodes!=masterVertexCount,std::runtime_error,
1754 "In call to STK_Interface::getElementVertices cardinality of " 1755 "element node relations must be the vertex count!");
1756 for(std::size_t node = 0; node < num_nodes; ++node) {
1760 for(
unsigned d=0;d<dim;d++)
1761 vertices_h(cell,node,d) = coord[d];
1764 Kokkos::deep_copy(vertices, vertices_h);
1767 template <
typename ArrayT>
1771 if(elements.size()==0) {
1780 unsigned masterVertexCount
1781 = stk::mesh::get_cell_topology(
bulkData_->bucket(elements[0]).topology()).getCellTopologyData()->vertex_count;
1785 auto vertices_h = Kokkos::create_mirror_view(vertices);
1786 for(std::size_t cell=0;cell<elements.size();cell++) {
1787 stk::mesh::Entity element = elements[cell];
1788 TEUCHOS_ASSERT(element!=0);
1790 unsigned vertexCount
1791 = stk::mesh::get_cell_topology(
bulkData_->bucket(element).topology()).getCellTopologyData()->vertex_count;
1792 TEUCHOS_TEST_FOR_EXCEPTION(vertexCount!=masterVertexCount,std::runtime_error,
1793 "In call to STK_Interface::getElementVertices all elements " 1794 "must have the same vertex count!");
1797 const size_t num_nodes =
bulkData_->num_nodes(element);
1798 stk::mesh::Entity
const* nodes =
bulkData_->begin_nodes(element);
1799 TEUCHOS_TEST_FOR_EXCEPTION(num_nodes!=masterVertexCount,std::runtime_error,
1800 "In call to STK_Interface::getElementVertices cardinality of " 1801 "element node relations must be the vertex count!");
1802 for(std::size_t node=0; node<num_nodes; ++node) {
1806 for(
unsigned d=0;d<dim;d++)
1807 vertices_h(cell,node,d) = coord[d];
1810 Kokkos::deep_copy(vertices, vertices_h);
1813 template <
typename ArrayT>
1819 if(elements.size()==0) {
1825 unsigned masterVertexCount
1826 = stk::mesh::get_cell_topology(
bulkData_->bucket(elements[0]).topology()).getCellTopologyData()->vertex_count;
1830 auto vertices_h = Kokkos::create_mirror_view(vertices);
1831 std::map<std::string,std::vector<std::string> >::const_iterator itr =
meshCoordFields_.find(eBlock);
1834 TEUCHOS_ASSERT(
false);
1837 const std::vector<std::string> & coordField = itr->second;
1838 std::vector<SolutionFieldType*> fields(
getDimension());
1839 for(std::size_t d=0;d<fields.size();d++) {
1843 for(std::size_t cell=0;cell<elements.size();cell++) {
1844 stk::mesh::Entity element = elements[cell];
1847 const size_t num_nodes =
bulkData_->num_nodes(element);
1848 stk::mesh::Entity
const* nodes =
bulkData_->begin_nodes(element);
1849 for(std::size_t i=0; i<num_nodes; ++i) {
1850 stk::mesh::Entity node = nodes[i];
1855 double * solnData = stk::mesh::field_data(*fields[d],node);
1859 vertices_h(cell,i,d) = solnData[0]+coord[d];
1863 Kokkos::deep_copy(vertices, vertices_h);
1866 template <
typename ArrayT>
1868 const std::string & eBlock, ArrayT & vertices)
const 1873 if(elements.size()==0) {
1877 std::map<std::string,std::vector<std::string> >::const_iterator itr =
meshCoordFields_.find(eBlock);
1880 TEUCHOS_ASSERT(
false);
1883 const std::vector<std::string> & coordField = itr->second;
1884 std::vector<SolutionFieldType*> fields(
getDimension());
1885 for(std::size_t d=0;d<fields.size();d++) {
1889 for(std::size_t cell=0;cell<elements.size();cell++) {
1890 stk::mesh::Entity element = elements[cell];
1893 const size_t num_nodes =
bulkData_->num_nodes(element);
1894 stk::mesh::Entity
const* nodes =
bulkData_->begin_nodes(element);
1895 for(std::size_t i=0; i<num_nodes; ++i) {
1896 stk::mesh::Entity node = nodes[i];
1901 double * solnData = stk::mesh::field_data(*fields[d],node);
1905 vertices(cell,i,d) = solnData[0]+coord[d];
std::size_t edgeLocalId(stk::mesh::Entity elmt) const
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 part
Teuchos::RCP< std::vector< stk::mesh::Entity > > orderedEdgeVector_
stk::mesh::EntityId faceGlobalId(std::size_t lid) const
bool isMeshCoordField(const std::string &eBlock, const std::string &fieldName, int &axis) const
ElementBlockException(const std::string &what)
void addNodeset(const std::string &name)
EdgeBlockException(const std::string &what)
void getElementVertices(const std::vector< std::size_t > &localIds, ArrayT &vertices) const
Teuchos::RCP< stk::mesh::BulkData > bulkData_
void setEdgeFieldData(const std::string &fieldName, const std::string &blockId, const std::vector< std::size_t > &localEdgeIds, const ArrayT &edgeValues, double scaleValue=1.0)
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)
bool isFaceLocal(stk::mesh::Entity face) const
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, Teuchos::RCP< shards::CellTopology > > faceBlockCT_
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)
std::map< std::string, stk::mesh::Part * > faceBlocks_
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_
void addEntityToFaceBlock(stk::mesh::Entity entity, stk::mesh::Part *faceblock)
stk::mesh::Field< double > * getEdgeField(const std::string &fieldName, const std::string &blockId) const
Teuchos::RCP< stk::mesh::BulkData > getBulkData() const
void getFaceBlockNames(std::vector< std::string > &names) const
std::size_t currentLocalId_
bool getUseFieldCoordinates() const
std::map< std::pair< std::string, std::string >, SolutionFieldType * > fieldNameToFaceField_
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
void addInformationRecords(const std::vector< std::string > &info_records)
void getMyFaces(std::vector< stk::mesh::Entity > &faces) 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)
stk::mesh::Part * getEdgeBlock(const std::string &name) const
get the block part
void getElementBlockNames(std::vector< std::string > &names) const
std::set< std::string > informationRecords_
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.
void setupExodusFile(const std::string &filename, const bool append=false, const bool append_after_restart_time=false, const double restart_time=0.0)
Set up an output Exodus file for writing results.
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)
FaceBlockException(const std::string &what)
static const std::string coordsString
stk::mesh::Field< double, stk::mesh::Cartesian > VectorFieldType
void setFaceFieldData(const std::string &fieldName, const std::string &blockId, const std::vector< std::size_t > &localFaceIds, const ArrayT &faceValues, double scaleValue=1.0)
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
std::map< std::string, stk::mesh::Part * > edgeBlocks_
void getElementVertices_FromCoords(const std::vector< stk::mesh::Entity > &elements, ArrayT &vertices) const
stk::mesh::Part * getFaceBlock(const std::string &name) const
get the block part
stk::mesh::EntityId edgeGlobalId(std::size_t lid) const
void getElementVertices_FromFieldNoResize(const std::vector< stk::mesh::Entity > &elements, const std::string &eBlock, ArrayT &vertices) const
Teuchos::RCP< const std::vector< stk::mesh::Entity > > getFacesOrderedByLID() const
std::map< std::string, stk::mesh::Part * > elementBlocks_
bool isInitialized() const
Has initialize been called on this mesh object?
stk::mesh::EntityId edgeGlobalId(stk::mesh::Entity edge) const
std::vector< stk::mesh::Part * > nodesPartVec_
static const std::string edgeBlockString
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_
void addFaceBlock(const std::string &name, const CellTopologyData *ctData)
Teuchos::RCP< const std::vector< stk::mesh::Entity > > getElementsOrderedByLID() const
void addNode(stk::mesh::EntityId gid, const std::vector< double > &coord)
void getAllFaces(const std::string &faceBlockName, std::vector< stk::mesh::Entity > &faces) const
void addEntitiesToEdgeBlock(std::vector< stk::mesh::Entity > entities, stk::mesh::Part *edgeblock)
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 getAllEdges(const std::string &edgeBlockName, std::vector< stk::mesh::Entity > &edges) const
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_
void getEdgeBlockNames(std::vector< std::string > &names) const
std::map< std::string, Teuchos::RCP< shards::CellTopology > > edgeBlockCT_
Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
get the comm associated with this mesh
void writeToExodus(const std::string &filename, const bool append=false)
Write this mesh and associated fields to the given output file.
void addEdgeField(const std::string &fieldName, const std::string &blockId)
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
static const std::string faceBlockString
std::unordered_map< stk::mesh::EntityId, std::size_t > localFaceIDHash_
void addEntityToEdgeBlock(stk::mesh::Entity entity, stk::mesh::Part *edgeblock)
void addEdgeBlock(const std::string &name, const CellTopologyData *ctData)
std::map< std::pair< std::string, std::string >, SolutionFieldType * > fieldNameToEdgeField_
void setCellFieldData(const std::string &fieldName, const std::string &blockId, const std::vector< std::size_t > &localElementIds, const ArrayT &solutionValues, double scaleValue=1.0)
std::unordered_map< stk::mesh::EntityId, std::size_t > localEdgeIDHash_
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
Teuchos::RCP< const std::vector< stk::mesh::Entity > > getEdgesOrderedByLID() 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
void addEntitiesToFaceBlock(std::vector< stk::mesh::Entity > entities, stk::mesh::Part *faceblock)
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
bool isEdgeLocal(stk::mesh::Entity edge) const
void refineMesh(const int numberOfLevels, const bool deleteParentElements)
void getNodesetNames(std::vector< std::string > &name) const
void buildLocalElementIDs()
Teuchos::RCP< std::vector< stk::mesh::Entity > > orderedFaceVector_
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
void getMyEdges(std::vector< stk::mesh::Entity > &edges) 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 addFaceField(const std::string &fieldName, const std::string &blockId)
void addCellField(const std::string &fieldName, const std::string &blockId)
stk::mesh::EntityId faceGlobalId(stk::mesh::Entity face) const
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
std::size_t faceLocalId(stk::mesh::Entity elmt) const
double getCurrentStateTime() const
stk::mesh::Part * getNodeset(const std::string &name) const
stk::mesh::Field< double > * getFaceField(const std::string &fieldName, const std::string &blockId) 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