43 #ifndef __Panzer_STK_Interface_hpp__
44 #define __Panzer_STK_Interface_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>
55 #include "Kokkos_Core.hpp"
57 #include <Shards_CellTopology.hpp>
58 #include <Shards_CellTopologyData.h>
60 #include <PanzerAdaptersSTK_config.hpp>
61 #include <Kokkos_ViewFactory.hpp>
63 #include <unordered_map>
65 #ifdef PANZER_HAVE_IOSS
66 #include <stk_io/StkMeshIoBroker.hpp>
69 #ifdef PANZER_HAVE_PERCEPT
72 class URP_Heterogeneous_3D;
76 namespace panzer_stk {
78 class PeriodicBC_MatcherBase;
85 ElementDescriptor(stk::mesh::EntityId gid,
const std::vector<stk::mesh::EntityId> & nodes);
92 std::vector<stk::mesh::EntityId>
nodes_;
135 void addElementBlock(
const std::string & name,
const CellTopologyData * ctData);
141 const std::string & edgeBlockName,
142 const stk::topology & topology);
144 const std::string & edgeBlockName,
145 const CellTopologyData * ctData);
150 const std::string & faceBlockName,
151 const stk::topology & topology);
153 const std::string & faceBlockName,
154 const CellTopologyData * ctData);
158 void addSideset(
const std::string & name,
const CellTopologyData * ctData);
166 void addSolutionField(
const std::string & fieldName,
const std::string & blockId);
170 void addCellField(
const std::string & fieldName,
const std::string & blockId);
174 void addEdgeField(
const std::string & fieldName,
const std::string & blockId);
178 void addFaceField(
const std::string & fieldName,
const std::string & blockId);
189 const std::vector<std::string> & coordField,
190 const std::string & dispPrefix);
211 void initialize(stk::ParallelMachine parallelMach,
bool setupIO=
true,
212 const bool buildRefinementSupport =
false);
237 void addNode(stk::mesh::EntityId gid,
const std::vector<double> & coord);
294 std::vector<stk::mesh::EntityId> & subcellIds)
const;
298 void getMyElements(std::vector<stk::mesh::Entity> & elements)
const;
302 void getMyElements(
const std::string & blockID,std::vector<stk::mesh::Entity> & elements)
const;
311 void getNeighborElements(
const std::string & blockID,std::vector<stk::mesh::Entity> & elements)
const;
315 void getMyEdges(std::vector<stk::mesh::Entity> & edges)
const;
324 void getMyEdges(
const std::string & edgeBlockName,std::vector<stk::mesh::Entity> & edges)
const;
334 void getMyEdges(
const std::string & edgeBlockName,
const std::string & blockName,std::vector<stk::mesh::Entity> & edges)
const;
343 void getAllEdges(
const std::string & edgeBlockName,std::vector<stk::mesh::Entity> & edges)
const;
353 void getAllEdges(
const std::string & edgeBlockName,
const std::string & blockName,std::vector<stk::mesh::Entity> & edges)
const;
357 void getMyFaces(std::vector<stk::mesh::Entity> & faces)
const;
366 void getMyFaces(
const std::string & faceBlockName,std::vector<stk::mesh::Entity> & faces)
const;
376 void getMyFaces(
const std::string & faceBlockName,
const std::string & blockName,std::vector<stk::mesh::Entity> & faces)
const;
385 void getAllFaces(
const std::string & faceBlockName,std::vector<stk::mesh::Entity> & faces)
const;
395 void getAllFaces(
const std::string & faceBlockName,
const std::string & blockName,std::vector<stk::mesh::Entity> & faces)
const;
404 void getMySides(
const std::string & sideName,std::vector<stk::mesh::Entity> & sides)
const;
414 void getMySides(
const std::string & sideName,
const std::string & blockName,std::vector<stk::mesh::Entity> & sides)
const;
423 void getAllSides(
const std::string & sideName,std::vector<stk::mesh::Entity> & sides)
const;
433 void getAllSides(
const std::string & sideName,
const std::string & blockName,std::vector<stk::mesh::Entity> & sides)
const;
444 void getMyNodes(
const std::string & sideName,
const std::string & blockName,std::vector<stk::mesh::Entity> & nodes)
const;
454 stk::mesh::Entity
findConnectivityById(stk::mesh::Entity src, stk::mesh::EntityRank tgt_rank,
unsigned rel_id)
const;
469 const bool append =
false);
490 const bool append =
false,
491 const bool append_after_restart_time =
false,
492 const double restart_time = 0.0);
531 const std::string& key,
555 const std::string& key,
556 const double& value);
579 const std::string& key,
580 const std::vector<int>& value);
603 const std::string& key,
604 const std::vector<double>& value);
615 #ifdef PANZER_HAVE_PERCEPT
624 {
if(
bulkData_==Teuchos::null)
return false;
625 return bulkData_->in_modifiable_state(); }
687 std::map<std::string, stk::mesh::Part*>::const_iterator itr =
elementBlocks_.find(name);
695 std::map<std::string, stk::mesh::Part*>::const_iterator itr =
edgeBlocks_.find(name);
703 std::map<std::string, stk::mesh::Part*>::const_iterator itr =
faceBlocks_.find(name);
715 return (itr !=
sidesets_.end()) ? itr->second :
nullptr;
725 return (itr !=
nodesets_.end()) ? itr->second :
nullptr;
741 void getNodeIdsForElement(stk::mesh::Entity element,std::vector<stk::mesh::EntityId> & nodeIds)
const;
746 std::vector<int> & relIds)
const;
751 std::vector<int> & relIds,
unsigned int matchType)
const;
755 void getElementsSharingNodes(
const std::vector<stk::mesh::EntityId> nodeId,std::vector<stk::mesh::Entity> & elements)
const;
788 std::size_t
edgeLocalId(stk::mesh::Entity elmt)
const;
792 std::size_t
edgeLocalId(stk::mesh::EntityId gid)
const;
814 std::size_t
faceLocalId(stk::mesh::Entity elmt)
const;
818 std::size_t
faceLocalId(stk::mesh::EntityId gid)
const;
833 {
return bulkData_->parallel_owner_rank(entity); }
837 inline bool isValid(stk::mesh::Entity entity)
const
849 const std::string & blockId)
const;
855 stk::mesh::Field<double> *
getCellField(
const std::string & fieldName,
856 const std::string & blockId)
const;
862 stk::mesh::Field<double> *
getEdgeField(
const std::string & fieldName,
863 const std::string & blockId)
const;
869 stk::mesh::Field<double> *
getFaceField(
const std::string & fieldName,
870 const std::string & blockId)
const;
896 template <
typename ArrayT>
898 const std::vector<std::size_t> & localElementIds,
const ArrayT & solutionValues,
double scaleValue=1.0);
914 template <
typename ArrayT>
916 const std::vector<std::size_t> & localElementIds,ArrayT & solutionValues)
const;
932 template <
typename ArrayT>
933 void setCellFieldData(
const std::string & fieldName,
const std::string & blockId,
934 const std::vector<std::size_t> & localElementIds,
const ArrayT & solutionValues,
double scaleValue=1.0);
960 template <
typename ArrayT>
961 void setEdgeFieldData(
const std::string & fieldName,
const std::string & blockId,
962 const std::vector<std::size_t> & localEdgeIds,
const ArrayT & edgeValues,
double scaleValue=1.0);
978 template <
typename ArrayT>
979 void setFaceFieldData(
const std::string & fieldName,
const std::string & blockId,
980 const std::vector<std::size_t> & localFaceIds,
const ArrayT & faceValues,
double scaleValue=1.0);
992 template <
typename ArrayT>
993 void getElementVertices(
const std::vector<std::size_t> & localIds, ArrayT & vertices)
const;
1003 template <
typename ArrayT>
1004 void getElementVertices(
const std::vector<stk::mesh::Entity> & elements, ArrayT & vertices)
const;
1015 template <
typename ArrayT>
1016 void getElementVertices(
const std::vector<std::size_t> & localIds,
const std::string & eBlock, ArrayT & vertices)
const;
1027 template <
typename ArrayT>
1028 void getElementVertices(
const std::vector<stk::mesh::Entity> & elements,
const std::string & eBlock, ArrayT & vertices)
const;
1038 template <
typename ArrayT>
1049 template <
typename ArrayT>
1061 template <
typename ArrayT>
1062 void getElementVerticesNoResize(
const std::vector<std::size_t> & localIds,
const std::string & eBlock, ArrayT & vertices)
const;
1073 template <
typename ArrayT>
1074 void getElementVerticesNoResize(
const std::vector<stk::mesh::Entity> & elements,
const std::string & eBlock, ArrayT & vertices)
const;
1087 template <
typename ArrayT>
1088 void getElementNodes(
const std::vector<std::size_t> & localIds, ArrayT & nodes)
const;
1098 template <
typename ArrayT>
1099 void getElementNodes(
const std::vector<stk::mesh::Entity> & elements, ArrayT & nodes)
const;
1110 template <
typename ArrayT>
1111 void getElementNodes(
const std::vector<std::size_t> & localIds,
const std::string & eBlock, ArrayT & nodes)
const;
1122 template <
typename ArrayT>
1123 void getElementNodes(
const std::vector<stk::mesh::Entity> & elements,
const std::string & eBlock, ArrayT & nodes)
const;
1133 template <
typename ArrayT>
1144 template <
typename ArrayT>
1156 template <
typename ArrayT>
1157 void getElementNodesNoResize(
const std::vector<std::size_t> & localIds,
const std::string & eBlock, ArrayT & nodes)
const;
1168 template <
typename ArrayT>
1169 void getElementNodesNoResize(
const std::vector<stk::mesh::Entity> & elements,
const std::string & eBlock, ArrayT & nodes)
const;
1176 stk::mesh::EntityRank
getFaceRank()
const {
return stk::topology::FACE_RANK; }
1177 stk::mesh::EntityRank
getEdgeRank()
const {
return stk::topology::EDGE_RANK; }
1178 stk::mesh::EntityRank
getNodeRank()
const {
return stk::topology::NODE_RANK; }
1198 const std::vector<Teuchos::RCP<const PeriodicBC_MatcherBase> > &
1204 std::vector<Teuchos::RCP<const PeriodicBC_MatcherBase> > &
1244 void print(std::ostream & os)
const;
1316 template <
typename ArrayT>
1317 void getElementVertices_FromField(
const std::vector<stk::mesh::Entity> & elements,
const std::string & eBlock, ArrayT & vertices)
const;
1319 template <
typename ArrayT>
1321 const std::string & eBlock, ArrayT & vertices)
const;
1332 template <
typename ArrayT>
1335 template <
typename ArrayT>
1350 template <
typename ArrayT>
1351 void getElementNodes_FromField(
const std::vector<stk::mesh::Entity> & elements,
const std::string & eBlock, ArrayT & nodes)
const;
1353 template <
typename ArrayT>
1355 const std::string & eBlock, ArrayT & nodes)
const;
1366 template <
typename ArrayT>
1369 template <
typename ArrayT>
1377 void refineMesh(
const int numberOfLevels,
const bool deleteParentElements);
1419 bool isMeshCoordField(
const std::string & eBlock,
const std::string & fieldName,
int & axis)
const;
1436 template <
typename ArrayT>
1437 void setDispFieldData(
const std::string & fieldName,
const std::string & blockId,
int axis,
1438 const std::vector<std::size_t> & localElementIds,
const ArrayT & solutionValues);
1445 #ifdef PANZER_HAVE_PERCEPT
1499 #ifdef PANZER_HAVE_IOSS
1508 enum class GlobalVariable
1532 const GlobalVariable& flag);
1581 template <
typename ArrayT>
1583 const std::vector<std::size_t> & localElementIds,
const ArrayT & solutionValues,
double scaleValue)
1586 auto solutionValues_h = Kokkos::create_mirror_view(solutionValues);
1587 Kokkos::deep_copy(solutionValues_h, solutionValues);
1589 int field_axis = -1;
1591 setDispFieldData(fieldName,blockId,field_axis,localElementIds,solutionValues_h);
1597 for(std::size_t cell=0;cell<localElementIds.size();cell++) {
1598 std::size_t localId = localElementIds[cell];
1599 stk::mesh::Entity element = elements[localId];
1602 const size_t num_nodes =
bulkData_->num_nodes(element);
1603 stk::mesh::Entity
const* nodes =
bulkData_->begin_nodes(element);
1604 for(std::size_t i=0; i<num_nodes; ++i) {
1605 stk::mesh::Entity node = nodes[i];
1607 double * solnData = stk::mesh::field_data(*field,node);
1609 solnData[0] = scaleValue*solutionValues_h(cell,i);
1614 template <
typename ArrayT>
1616 const std::vector<std::size_t> & localElementIds,
const ArrayT & dispValues)
1625 for(std::size_t cell=0;cell<localElementIds.size();cell++) {
1626 std::size_t localId = localElementIds[cell];
1627 stk::mesh::Entity element = elements[localId];
1630 const size_t num_nodes =
bulkData_->num_nodes(element);
1631 stk::mesh::Entity
const* nodes =
bulkData_->begin_nodes(element);
1632 for(std::size_t i=0; i<num_nodes; ++i) {
1633 stk::mesh::Entity node = nodes[i];
1635 double * solnData = stk::mesh::field_data(*field,node);
1636 double * coordData = stk::mesh::field_data(coord_field,node);
1638 solnData[0] = dispValues(cell,i)-coordData[axis];
1643 template <
typename ArrayT>
1645 const std::vector<std::size_t> & localElementIds,ArrayT & solutionValues)
const
1651 localElementIds.size(),
1652 bulkData_->num_nodes(elements[localElementIds[0]]));
1656 for(std::size_t cell=0;cell<localElementIds.size();cell++) {
1657 std::size_t localId = localElementIds[cell];
1658 stk::mesh::Entity element = elements[localId];
1661 const size_t num_nodes =
bulkData_->num_nodes(element);
1662 stk::mesh::Entity
const* nodes =
bulkData_->begin_nodes(element);
1663 for(std::size_t i=0; i<num_nodes; ++i) {
1664 stk::mesh::Entity node = nodes[i];
1666 double * solnData = stk::mesh::field_data(*field,node);
1668 solutionValues(cell,i) = solnData[0];
1673 template <
typename ArrayT>
1675 const std::vector<std::size_t> & localElementIds,
const ArrayT & solutionValues,
double scaleValue)
1681 auto solutionValues_h = Kokkos::create_mirror_view(solutionValues);
1682 Kokkos::deep_copy(solutionValues_h, solutionValues);
1684 for(std::size_t cell=0;cell<localElementIds.size();cell++) {
1685 std::size_t localId = localElementIds[cell];
1686 stk::mesh::Entity element = elements[localId];
1688 double * solnData = stk::mesh::field_data(*field,element);
1690 solnData[0] = scaleValue*solutionValues_h.access(cell,0);
1694 template <
typename ArrayT>
1696 const std::vector<std::size_t> & localEdgeIds,
const ArrayT & edgeValues,
double scaleValue)
1702 auto edgeValues_h = Kokkos::create_mirror_view(edgeValues);
1703 Kokkos::deep_copy(edgeValues_h, edgeValues);
1705 for(std::size_t idx=0;idx<localEdgeIds.size();idx++) {
1706 std::size_t localId = localEdgeIds[idx];
1707 stk::mesh::Entity edge = edges[localId];
1709 double * solnData = stk::mesh::field_data(*field,edge);
1711 solnData[0] = scaleValue*edgeValues_h.access(idx,0);
1715 template <
typename ArrayT>
1717 const std::vector<std::size_t> & localFaceIds,
const ArrayT & faceValues,
double scaleValue)
1723 auto faceValues_h = Kokkos::create_mirror_view(faceValues);
1724 Kokkos::deep_copy(faceValues_h, faceValues);
1726 for(std::size_t idx=0;idx<localFaceIds.size();idx++) {
1727 std::size_t localId = localFaceIds[idx];
1728 stk::mesh::Entity face = faces[localId];
1730 double * solnData = stk::mesh::field_data(*field,face);
1732 solnData[0] = scaleValue*faceValues_h.access(idx,0);
1738 template <
typename ArrayT>
1749 std::vector<stk::mesh::Entity> selected_elements;
1750 for(std::size_t cell=0;cell<localElementIds.size();cell++)
1751 selected_elements.push_back(elements[localElementIds[cell]]);
1757 "STK_Interface::getElementVertices: Cannot call this method when field coordinates are used "
1758 "without specifying an element block.");
1762 template <
typename ArrayT>
1770 "STK_Interface::getElementVertices: Cannot call this method when field coordinates are used "
1771 "without specifying an element block.");
1775 template <
typename ArrayT>
1786 template <
typename ArrayT>
1792 std::vector<stk::mesh::Entity> selected_elements;
1793 for(std::size_t cell=0;cell<localElementIds.size();cell++)
1794 selected_elements.push_back(elements[localElementIds[cell]]);
1804 template <
typename ArrayT>
1815 std::vector<stk::mesh::Entity> selected_elements;
1816 for(std::size_t cell=0;cell<localElementIds.size();cell++)
1817 selected_elements.push_back(elements[localElementIds[cell]]);
1823 "STK_Interface::getElementVerticesNoResize: Cannot call this method when field coordinates are used "
1824 "without specifying an element block.");
1828 template <
typename ArrayT>
1836 "STK_Interface::getElementVerticesNoResize: Cannot call this method when field coordinates are used "
1837 "without specifying an element block.");
1841 template <
typename ArrayT>
1852 template <
typename ArrayT>
1858 std::vector<stk::mesh::Entity> selected_elements;
1859 for(std::size_t cell=0;cell<localElementIds.size();cell++)
1860 selected_elements.push_back(elements[localElementIds[cell]]);
1870 template <
typename ArrayT>
1874 if(elements.size() == 0) {
1884 const auto masterVertexCount
1885 = stk::mesh::get_cell_topology(
bulkData_->bucket(elements[0]).topology()).getCellTopologyData()->vertex_count;
1889 auto vertices_h = Kokkos::create_mirror_view(vertices);
1890 Kokkos::deep_copy(vertices_h, vertices);
1894 for(std::size_t cell = 0; cell < elements.size(); cell++) {
1895 const auto element = elements[cell];
1898 const auto vertexCount
1899 = stk::mesh::get_cell_topology(
bulkData_->bucket(element).topology()).getCellTopologyData()->vertex_count;
1901 "In call to STK_Interface::getElementVertices all elements "
1902 "must have the same vertex count!");
1905 const size_t num_nodes =
bulkData_->num_nodes(element);
1906 auto const* nodes =
bulkData_->begin_nodes(element);
1908 "In call to STK_Interface::getElementVertices cardinality of "
1909 "element node relations must be the vertex count!");
1910 for(std::size_t node = 0; node < num_nodes; ++node) {
1914 for(
unsigned d=0;d<dim;d++)
1915 vertices_h(cell,node,d) = coord[d];
1918 Kokkos::deep_copy(vertices, vertices_h);
1921 template <
typename ArrayT>
1925 if(elements.size()==0) {
1934 unsigned masterVertexCount
1935 = stk::mesh::get_cell_topology(
bulkData_->bucket(elements[0]).topology()).getCellTopologyData()->vertex_count;
1939 auto vertices_h = Kokkos::create_mirror_view(vertices);
1940 for(std::size_t cell=0;cell<elements.size();cell++) {
1941 stk::mesh::Entity element = elements[cell];
1944 unsigned vertexCount
1945 = stk::mesh::get_cell_topology(
bulkData_->bucket(element).topology()).getCellTopologyData()->vertex_count;
1947 "In call to STK_Interface::getElementVertices all elements "
1948 "must have the same vertex count!");
1951 const size_t num_nodes =
bulkData_->num_nodes(element);
1952 stk::mesh::Entity
const* nodes =
bulkData_->begin_nodes(element);
1954 "In call to STK_Interface::getElementVertices cardinality of "
1955 "element node relations must be the vertex count!");
1956 for(std::size_t node=0; node<num_nodes; ++node) {
1960 for(
unsigned d=0;d<dim;d++)
1961 vertices_h(cell,node,d) = coord[d];
1964 Kokkos::deep_copy(vertices, vertices_h);
1967 template <
typename ArrayT>
1973 if(elements.size()==0) {
1979 unsigned masterVertexCount
1980 = stk::mesh::get_cell_topology(
bulkData_->bucket(elements[0]).topology()).getCellTopologyData()->vertex_count;
1984 auto vertices_h = Kokkos::create_mirror_view(vertices);
1985 std::map<std::string,std::vector<std::string> >::const_iterator itr =
meshCoordFields_.find(eBlock);
1991 const std::vector<std::string> & coordField = itr->second;
1992 std::vector<SolutionFieldType*> fields(
getDimension());
1993 for(std::size_t d=0;d<fields.size();d++) {
1997 for(std::size_t cell=0;cell<elements.size();cell++) {
1998 stk::mesh::Entity element = elements[cell];
2001 const size_t num_nodes =
bulkData_->num_nodes(element);
2002 stk::mesh::Entity
const* nodes =
bulkData_->begin_nodes(element);
2003 for(std::size_t i=0; i<num_nodes; ++i) {
2004 stk::mesh::Entity node = nodes[i];
2009 double * solnData = stk::mesh::field_data(*fields[d],node);
2013 vertices_h(cell,i,d) = solnData[0]+coord[d];
2017 Kokkos::deep_copy(vertices, vertices_h);
2020 template <
typename ArrayT>
2022 const std::string & eBlock, ArrayT & vertices)
const
2027 if(elements.size()==0) {
2031 std::map<std::string,std::vector<std::string> >::const_iterator itr =
meshCoordFields_.find(eBlock);
2037 const std::vector<std::string> & coordField = itr->second;
2038 std::vector<SolutionFieldType*> fields(
getDimension());
2039 for(std::size_t d=0;d<fields.size();d++) {
2043 for(std::size_t cell=0;cell<elements.size();cell++) {
2044 stk::mesh::Entity element = elements[cell];
2047 const size_t num_nodes =
bulkData_->num_nodes(element);
2048 stk::mesh::Entity
const* nodes =
bulkData_->begin_nodes(element);
2049 for(std::size_t i=0; i<num_nodes; ++i) {
2050 stk::mesh::Entity node = nodes[i];
2055 double * solnData = stk::mesh::field_data(*fields[d],node);
2059 vertices(cell,i,d) = solnData[0]+coord[d];
2067 template <
typename ArrayT>
2078 std::vector<stk::mesh::Entity> selected_elements;
2079 for(std::size_t cell=0;cell<localElementIds.size();cell++)
2080 selected_elements.push_back(elements[localElementIds[cell]]);
2086 "STK_Interface::getElementNodes: Cannot call this method when field coordinates are used "
2087 "without specifying an element block.");
2091 template <
typename ArrayT>
2099 "STK_Interface::getElementNodes: Cannot call this method when field coordinates are used "
2100 "without specifying an element block.");
2104 template <
typename ArrayT>
2115 template <
typename ArrayT>
2121 std::vector<stk::mesh::Entity> selected_elements;
2122 for(std::size_t cell=0;cell<localElementIds.size();cell++)
2123 selected_elements.push_back(elements[localElementIds[cell]]);
2133 template <
typename ArrayT>
2144 std::vector<stk::mesh::Entity> selected_elements;
2145 for(std::size_t cell=0;cell<localElementIds.size();cell++)
2146 selected_elements.push_back(elements[localElementIds[cell]]);
2152 "STK_Interface::getElementNodesNoResize: Cannot call this method when field coordinates are used "
2153 "without specifying an element block.");
2157 template <
typename ArrayT>
2165 "STK_Interface::getElementNodesNoResize: Cannot call this method when field coordinates are used "
2166 "without specifying an element block.");
2170 template <
typename ArrayT>
2181 template <
typename ArrayT>
2187 std::vector<stk::mesh::Entity> selected_elements;
2188 for(std::size_t cell=0;cell<localElementIds.size();cell++)
2189 selected_elements.push_back(elements[localElementIds[cell]]);
2199 template <
typename ArrayT>
2203 if(elements.size() == 0) {
2213 const auto masterNodeCount
2214 = stk::mesh::get_cell_topology(
bulkData_->bucket(elements[0]).topology()).getCellTopologyData()->node_count;
2218 auto nodes_h = Kokkos::create_mirror_view(nodes);
2219 Kokkos::deep_copy(nodes_h, nodes);
2223 for(std::size_t cell = 0; cell < elements.size(); cell++) {
2224 const auto element = elements[cell];
2227 const auto nodeCount
2228 = stk::mesh::get_cell_topology(
bulkData_->bucket(element).topology()).getCellTopologyData()->node_count;
2230 "In call to STK_Interface::getElementNodes all elements "
2231 "must have the same node count!");
2234 const size_t num_nodes =
bulkData_->num_nodes(element);
2235 auto const* elem_nodes =
bulkData_->begin_nodes(element);
2237 "In call to STK_Interface::getElementNodes cardinality of "
2238 "element node relations must be the node count!");
2239 for(std::size_t node = 0; node < num_nodes; ++node) {
2243 for(
unsigned d=0;d<dim;d++)
2244 nodes_h(cell,node,d) = coord[d];
2247 Kokkos::deep_copy(nodes, nodes_h);
2250 template <
typename ArrayT>
2254 if(elements.size()==0) {
2263 unsigned masterNodeCount
2264 = stk::mesh::get_cell_topology(
bulkData_->bucket(elements[0]).topology()).getCellTopologyData()->node_count;
2268 auto nodes_h = Kokkos::create_mirror_view(nodes);
2269 for(std::size_t cell=0;cell<elements.size();cell++) {
2270 stk::mesh::Entity element = elements[cell];
2274 = stk::mesh::get_cell_topology(
bulkData_->bucket(element).topology()).getCellTopologyData()->node_count;
2276 "In call to STK_Interface::getElementNodes all elements "
2277 "must have the same node count!");
2280 const size_t num_nodes =
bulkData_->num_nodes(element);
2281 stk::mesh::Entity
const* elem_nodes =
bulkData_->begin_nodes(element);
2283 "In call to STK_Interface::getElementNodes cardinality of "
2284 "element node relations must be the node count!");
2285 for(std::size_t node=0; node<num_nodes; ++node) {
2289 for(
unsigned d=0;d<dim;d++)
2290 nodes_h(cell,node,d) = coord[d];
2293 Kokkos::deep_copy(nodes, nodes_h);
2296 template <
typename ArrayT>
2302 if(elements.size()==0) {
2308 unsigned masterNodeCount
2309 = stk::mesh::get_cell_topology(
bulkData_->bucket(elements[0]).topology()).getCellTopologyData()->node_count;
2313 auto nodes_h = Kokkos::create_mirror_view(nodes);
2314 std::map<std::string,std::vector<std::string> >::const_iterator itr =
meshCoordFields_.find(eBlock);
2320 const std::vector<std::string> & coordField = itr->second;
2321 std::vector<SolutionFieldType*> fields(
getDimension());
2322 for(std::size_t d=0;d<fields.size();d++) {
2327 for(std::size_t cell=0;cell<elements.size();cell++) {
2328 stk::mesh::Entity element = elements[cell];
2331 const size_t num_nodes =
bulkData_->num_nodes(element);
2332 stk::mesh::Entity
const* elem_nodes =
bulkData_->begin_nodes(element);
2333 for(std::size_t i=0; i<num_nodes; ++i) {
2334 stk::mesh::Entity node = elem_nodes[i];
2339 double * solnData = stk::mesh::field_data(*fields[d],node);
2343 nodes_h(cell,i,d) = solnData[0]+coord[d];
2347 Kokkos::deep_copy(nodes, nodes_h);
2350 template <
typename ArrayT>
2352 const std::string & eBlock, ArrayT & nodes)
const
2357 if(elements.size()==0) {
2361 std::map<std::string,std::vector<std::string> >::const_iterator itr =
meshCoordFields_.find(eBlock);
2367 const std::vector<std::string> & coordField = itr->second;
2368 std::vector<SolutionFieldType*> fields(
getDimension());
2369 for(std::size_t d=0;d<fields.size();d++) {
2374 for(std::size_t cell=0;cell<elements.size();cell++) {
2375 stk::mesh::Entity element = elements[cell];
2378 const size_t num_nodes =
bulkData_->num_nodes(element);
2379 stk::mesh::Entity
const* elem_nodes =
bulkData_->begin_nodes(element);
2380 for(std::size_t i=0; i<num_nodes; ++i) {
2381 stk::mesh::Entity node = elem_nodes[i];
2386 double * solnData = stk::mesh::field_data(*fields[d],node);
2390 nodes(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< Teuchos::MpiComm< int > > getSafeCommunicator(stk::ParallelMachine parallelMach) const
void getSidesetNames(std::vector< std::string > &name) const
std::pair< Teuchos::RCP< std::vector< std::pair< std::size_t, std::size_t > > >, Teuchos::RCP< std::vector< unsigned int > > > getPeriodicNodePairing() const
const bool & useBoundingBoxSearch() const
Teuchos::RCP< std::vector< stk::mesh::Entity > > orderedEdgeVector_
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.
ElementBlockException(const std::string &what)
void addNodeset(const std::string &name)
EdgeBlockException(const std::string &what)
Teuchos::RCP< stk::mesh::BulkData > bulkData_
void getFaceBlockNames(std::vector< std::string > &names) const
void print(std::ostream &os) const
void getElementNodesNoResize(const std::vector< std::size_t > &localIds, ArrayT &nodes) const
unsigned entityOwnerRank(stk::mesh::Entity entity) const
void getMyEdges(std::vector< stk::mesh::Entity > &edges) const
void setEdgeFieldData(const std::string &fieldName, const std::string &blockId, const std::vector< std::size_t > &localEdgeIds, const ArrayT &edgeValues, double scaleValue=1.0)
stk::mesh::Field< double > VectorFieldType
Teuchos::RCP< const std::vector< stk::mesh::Entity > > getElementsOrderedByLID() const
void setSolutionFieldData(const std::string &fieldName, const std::string &blockId, const std::vector< std::size_t > &localElementIds, const ArrayT &solutionValues, double scaleValue=1.0)
std::vector< stk::mesh::Part * > facesPartVec_
std::vector< stk::mesh::EntityId > nodes_
std::map< std::pair< std::string, std::string >, SolutionFieldType * > fieldNameToSolution_
void getElementVerticesNoResize(const std::vector< std::size_t > &localIds, ArrayT &vertices) const
void addEntityToNodeset(stk::mesh::Entity entity, stk::mesh::Part *nodeset)
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 getElementNodes_FromCoordsNoResize(const std::vector< stk::mesh::Entity > &elements, ArrayT &nodes) const
void getElementBlockNames(std::vector< std::string > &names) const
void addPeriodicBC(const Teuchos::RCP< const PeriodicBC_MatcherBase > &bc)
void getSubcellIndices(unsigned entityRank, stk::mesh::EntityId elementId, std::vector< stk::mesh::EntityId > &subcellIds) const
std::vector< Teuchos::RCP< const PeriodicBC_MatcherBase > > & getPeriodicBCVector()
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::Field< double > SolutionFieldType
stk::mesh::Part * getElementBlockPart(const std::string &name) const
get the block part
std::vector< stk::mesh::Part * > edgesPartVec_
bool nonnull(const std::shared_ptr< T > &p)
bool isFaceLocal(stk::mesh::Entity face) const
std::map< std::string, double > blockWeights_
Teuchos::RCP< const std::vector< stk::mesh::Entity > > getEdgesOrderedByLID() const
VectorFieldType * coordinatesField_
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
std::size_t getEntityCounts(unsigned entityRank) const
get the global counts for the entity of specified rank
void printMetaData(std::ostream &os) const
void getElementNodes_FromFieldNoResize(const std::vector< stk::mesh::Entity > &elements, const std::string &eBlock, ArrayT &nodes) const
void getElementVertices(const std::vector< std::size_t > &localIds, ArrayT &vertices) const
void addEntityToSideset(stk::mesh::Entity entity, stk::mesh::Part *sideset)
std::map< std::string, stk::mesh::Part * > faceBlocks_
stk::mesh::Part * nodesPart_
stk::mesh::EntityId edgeGlobalId(std::size_t lid) const
void getElementVertices_FromCoords(const std::vector< stk::mesh::Entity > &elements, ArrayT &vertices) const
void getElementNodes(const std::vector< std::size_t > &localIds, ArrayT &nodes) 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::EntityId getMaxEntityId(unsigned entityRank) const
get max entity ID of type entityRank
std::size_t currentLocalId_
std::map< std::pair< std::string, std::string >, SolutionFieldType * > fieldNameToFaceField_
std::map< std::string, stk::mesh::Part * > sidesets_
std::vector< std::size_t > entityCounts_
bool validBlockId(const std::string &blockId) const
stk::mesh::Part * getNodeset(const std::string &name) const
void getEdgeBlockNames(std::vector< std::string > &names) const
stk::mesh::Part * getFaceBlock(const std::string &name) const
get the block part
void addInformationRecords(const std::vector< std::string > &info_records)
stk::mesh::EntityId faceGlobalId(std::size_t lid) const
void getAllFaces(const std::string &faceBlockName, std::vector< stk::mesh::Entity > &faces) const
stk::mesh::Field< double > * getSolutionField(const std::string &fieldName, const std::string &blockId) const
SolutionFieldType * loadBalField_
void addEdgeBlock(const std::string &elemBlockName, const std::string &edgeBlockName, const stk::topology &topology)
const double * getNodeCoordinates(stk::mesh::EntityId nodeId) const
stk::mesh::Part * getSideset(const std::string &name) const
static const std::string nodesString
void instantiateBulkData(stk::ParallelMachine parallelMach)
void getElementVertices_FromFieldNoResize(const std::vector< stk::mesh::Entity > &elements, const std::string &eBlock, ArrayT &vertices) const
std::set< std::string > informationRecords_
void getMySides(const std::string &sideName, std::vector< stk::mesh::Entity > &sides) const
std::size_t getNumSidesets() const
get the side set count
stk::mesh::EntityId elementGlobalId(stk::mesh::Entity elmt) const
stk::mesh::Field< double > * getEdgeField(const std::string &fieldName, const std::string &blockId) const
std::size_t elementLocalId(stk::mesh::Entity elmt) const
bool useFieldCoordinates_
stk::mesh::EntityId faceGlobalId(stk::mesh::Entity face) const
void getElementVertices_FromCoordsNoResize(const std::vector< stk::mesh::Entity > &elements, ArrayT &vertices) const
unsigned getDimension() const
get the dimension
void addElement(const Teuchos::RCP< ElementDescriptor > &ed, stk::mesh::Part *block)
std::map< std::string, std::vector< std::string > > meshDispFields_
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::EntityRank getSideRank() const
void getAllEdges(const std::string &edgeBlockName, std::vector< stk::mesh::Entity > &edges) const
void initialize(stk::ParallelMachine parallelMach, bool setupIO=true, const bool buildRefinementSupport=false)
FaceBlockException(const std::string &what)
static const std::string coordsString
void setFaceFieldData(const std::string &fieldName, const std::string &blockId, const std::vector< std::size_t > &localFaceIds, const ArrayT &faceValues, double scaleValue=1.0)
std::size_t getNumElementBlocks() const
get the block count
void getElementNodes_FromField(const std::vector< stk::mesh::Entity > &elements, const std::string &eBlock, ArrayT &nodes) const
Teuchos::RCP< stk::mesh::BulkData > getBulkData() const
std::map< std::string, stk::mesh::Part * > edgeBlocks_
void getElementVertices_FromField(const std::vector< stk::mesh::Entity > &elements, const std::string &eBlock, ArrayT &vertices) const
bool getUseFieldCoordinates() const
stk::mesh::EntityId elementGlobalId(std::size_t lid) const
std::map< std::string, stk::mesh::Part * > elementBlocks_
std::vector< stk::mesh::Part * > nodesPartVec_
static const std::string edgeBlockString
static const std::string edgesString
SidesetException(const std::string &what)
std::string containingBlockId(stk::mesh::Entity elmt) const
stk::mesh::EntityRank getFaceRank() const
std::size_t edgeLocalId(stk::mesh::Entity elmt) const
const STK_Interface * mesh_
std::vector< Teuchos::RCP< const PeriodicBC_MatcherBase > > periodicBCs_
void getOwnedElementsSharingNode(stk::mesh::Entity node, std::vector< stk::mesh::Entity > &elements, std::vector< int > &relIds) const
void addNode(stk::mesh::EntityId gid, const std::vector< double > &coord)
void addEntitiesToEdgeBlock(std::vector< stk::mesh::Entity > entities, stk::mesh::Part *edgeblock)
void setBlockWeight(const std::string &blockId, double weight)
void buildSubcells()
force the mesh to build subcells: edges and faces
Teuchos::RCP< Teuchos::MpiComm< int > > mpiComm_
LocalIdCompare(const STK_Interface *mesh)
void setInitialStateTime(double value)
Teuchos::RCP< stk::mesh::MetaData > metaData_
Teuchos::RCP< stk::mesh::MetaData > getMetaData() const
virtual ~ElementDescriptor()
Teuchos::RCP< std::vector< stk::mesh::Entity > > orderedElementVector_
std::size_t getNumNodesets() const
get the side set count
double getInitialStateTime() const
bool isInitialized() const
Has initialize been called on this mesh object?
void writeToExodus(const std::string &filename, const bool append=false)
Write this mesh and associated fields to the given output file.
void getMyNodes(const std::string &sideName, const std::string &blockName, std::vector< stk::mesh::Entity > &nodes) const
void getSolutionFieldData(const std::string &fieldName, const std::string &blockId, const std::vector< std::size_t > &localElementIds, ArrayT &solutionValues) const
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()
stk::mesh::Entity findConnectivityById(stk::mesh::Entity src, stk::mesh::EntityRank tgt_rank, unsigned rel_id) const
void setUseLowerCaseForIO(bool useLowerCase)
bool isMeshCoordField(const std::string &eBlock, const std::string &fieldName, int &axis) const
void getAllSides(const std::string &sideName, std::vector< stk::mesh::Entity > &sides) const
void addSideset(const std::string &name, const CellTopologyData *ctData)
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)
std::map< std::pair< std::string, std::string >, SolutionFieldType * > fieldNameToEdgeField_
stk::mesh::Part * getEdgeBlock(const std::string &name) const
get the block part
void setCellFieldData(const std::string &fieldName, const std::string &blockId, const std::vector< std::size_t > &localElementIds, const ArrayT &solutionValues, double scaleValue=1.0)
const VectorFieldType & getEdgesField() const
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)
std::map< std::string, std::vector< std::string > > meshCoordFields_
stk::mesh::EntityRank getEdgeRank() const
static const std::string facesString
void addEntitiesToFaceBlock(std::vector< stk::mesh::Entity > entities, stk::mesh::Part *faceblock)
void initializeFieldsInSTK(const std::map< std::pair< std::string, std::string >, SolutionFieldType * > &nameToField, bool setupIO)
Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
get the comm associated with this mesh
const std::vector< Teuchos::RCP< const PeriodicBC_MatcherBase > > & getPeriodicBCVector() const
const VectorFieldType & getFacesField() const
void refineMesh(const int numberOfLevels, const bool deleteParentElements)
stk::mesh::Field< double > * getFaceField(const std::string &fieldName, const std::string &blockId) const
void buildLocalElementIDs()
std::size_t faceLocalId(stk::mesh::Entity elmt) const
bool isModifiable() const
void addFaceBlock(const std::string &elemBlockName, const std::string &faceBlockName, const stk::topology &topology)
Teuchos::RCP< std::vector< stk::mesh::Entity > > orderedFaceVector_
std::unordered_map< stk::mesh::EntityId, std::size_t > localIDHash_
bool isValid(stk::mesh::Entity entity) const
std::vector< stk::mesh::EntityId > maxEntityId_
stk::mesh::Field< ProcIdData > ProcIdFieldType
void addPeriodicBCs(const std::vector< Teuchos::RCP< const PeriodicBC_MatcherBase > > &bc_vec)
void getElementsSharingNode(stk::mesh::EntityId nodeId, std::vector< stk::mesh::Entity > &elements) const
get a set of elements sharing a single node
void addFaceField(const std::string &fieldName, const std::string &blockId)
#define TEUCHOS_ASSERT(assertion_test)
void getMyElements(std::vector< stk::mesh::Entity > &elements) const
stk::mesh::EntityRank getNodeRank() const
void getMyFaces(std::vector< stk::mesh::Entity > &faces) const
void addCellField(const std::string &fieldName, const std::string &blockId)
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
Teuchos::RCP< const std::vector< stk::mesh::Entity > > getFacesOrderedByLID() const
stk::mesh::Part * edgesPart_
bool operator()(stk::mesh::Entity a, stk::mesh::Entity b)
stk::mesh::EntityRank getElementRank() const
void getElementNodes_FromCoords(const std::vector< stk::mesh::Entity > &elements, ArrayT &nodes) const
stk::mesh::EntityId edgeGlobalId(stk::mesh::Entity edge) const
stk::mesh::Part * facesPart_
ProcIdFieldType * processorIdField_
void addElementBlock(const std::string &name, const CellTopologyData *ctData)
const std::vector< stk::mesh::EntityId > & getNodes() const
stk::mesh::Field< double > * getCellField(const std::string &fieldName, const std::string &blockId) const
const VectorFieldType & getCoordinatesField() const
stk::mesh::Part * getOwnedPart() const
Get a pointer to the locally owned part.
stk::mesh::EntityId getGID() const
void setBoundingBoxSearchFlag(const bool &searchFlag)
bool getUseLowerCaseForIO() const
bool isEdgeLocal(stk::mesh::Entity edge) const
double getCurrentStateTime() const
void getNeighborElements(std::vector< stk::mesh::Entity > &elements) const
void getNodesetNames(std::vector< std::string > &name) const
Teuchos::RCP< const shards::CellTopology > getCellTopology(const std::string &eBlock) const
std::map< std::pair< std::string, std::string >, SolutionFieldType * > fieldNameToCellField_
std::map< std::string, stk::mesh::Part * > nodesets_
VectorFieldType * facesField_