2 #include <stk_mesh/fem/FEMMetaData.hpp> 4 #include <Shards_CellTopology.hpp> 5 #include <stk_mesh/base/BulkData.hpp> 6 #include <stk_mesh/base/Bucket.hpp> 7 #include <stk_mesh/base/Entity.hpp> 8 #include <stk_mesh/base/Ghosting.hpp> 16 void assign_cell_topology(
19 const fem::CellTopology cell_topology)
21 if (part_ordinal >= part_cell_topology_vector.size())
22 part_cell_topology_vector.resize(part_ordinal + 1);
24 part_cell_topology_vector[part_ordinal] = cell_topology;
26 if (!cell_topology.getCellTopologyData())
28 std::cout <<
"bad topology in FEMMetaData::assign_cell_topology" << std::endl;
31 ThrowRequireMsg(cell_topology.getCellTopologyData(),
"bad topology in FEMMetaData::assign_cell_topology");
36 FEMMetaData::FEMMetaData()
38 m_fem_initialized(false),
39 m_spatial_dimension(0),
40 m_side_rank(INVALID_RANK),
41 m_element_rank(INVALID_RANK)
44 m_meta_data.declare_attribute_no_delete<FEMMetaData>(
this);
47 FEMMetaData::FEMMetaData(
size_t spatial_dimension,
48 const std::vector<std::string>& in_entity_rank_names)
50 m_fem_initialized(false),
51 m_spatial_dimension(0),
52 m_side_rank(INVALID_RANK),
53 m_element_rank(INVALID_RANK)
56 m_meta_data.declare_attribute_no_delete<
FEMMetaData>(
this);
63 ThrowRequireMsg(!m_fem_initialized,
"FEM functionality in FEMMetaData can only be initialized once.");
64 if ( rank_names.empty() ) {
69 "Entity rank name vector must name every rank");
70 m_entity_rank_names = rank_names;
74 m_fem_initialized =
true;
75 internal_declare_known_cell_topology_parts();
78 void FEMMetaData::internal_set_spatial_dimension_and_ranks(
size_t spatial_dimension)
80 ThrowRequireMsg(
spatial_dimension != 0,
"FEMMetaData::internal_set_spatial_dimension_and_ranks: spatial_dimension == 0!");
93 m_side_rank = m_spatial_dimension - 1;
94 m_element_rank = m_spatial_dimension;
98 void FEMMetaData::internal_declare_known_cell_topology_parts()
103 if (m_spatial_dimension == 1) {
105 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Particle >()), m_element_rank);
107 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Line<2> >()), m_element_rank);
108 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Line<3> >()), m_element_rank);
112 else if (m_spatial_dimension == 2) {
117 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Particle >()), m_element_rank);
119 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Triangle<3> >()), m_element_rank);
120 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Triangle<6> >()), m_element_rank);
121 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Triangle<4> >()), m_element_rank);
123 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Quadrilateral<4> >()), m_element_rank);
124 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Quadrilateral<8> >()), m_element_rank);
125 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Quadrilateral<9> >()), m_element_rank);
127 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Beam<2> >()), m_element_rank);
128 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Beam<3> >()), m_element_rank);
130 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::ShellLine<2> >()), m_element_rank);
131 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::ShellLine<3> >()), m_element_rank);
134 else if (m_spatial_dimension == 3) {
139 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Triangle<3> >()), m_side_rank);
140 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Triangle<6> >()), m_side_rank);
141 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Triangle<4> >()), m_side_rank);
143 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Quadrilateral<4> >()), m_side_rank);
144 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Quadrilateral<8> >()), m_side_rank);
145 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Quadrilateral<9> >()), m_side_rank);
147 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Particle >()), m_element_rank);
149 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Beam<2> >()), m_element_rank);
150 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Beam<3> >()), m_element_rank);
152 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Tetrahedron<4> >()), m_element_rank);
153 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Tetrahedron<10> >()), m_element_rank);
154 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Tetrahedron<11> >()), m_element_rank);
155 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Tetrahedron<8> >()), m_element_rank);
157 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Pyramid<5> >()), m_element_rank);
158 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Pyramid<13> >()), m_element_rank);
159 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Pyramid<14> >()), m_element_rank);
161 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Wedge<6> >()), m_element_rank);
162 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Wedge<15> >()), m_element_rank);
163 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Wedge<18> >()), m_element_rank);
165 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Hexahedron<8> >()), m_element_rank);
166 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Hexahedron<20> >()), m_element_rank);
167 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::Hexahedron<27> >()), m_element_rank);
169 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::ShellTriangle<3> >()), m_element_rank);
170 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::ShellTriangle<6> >()), m_element_rank);
172 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::ShellQuadrilateral<4> >()), m_element_rank);
173 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::ShellQuadrilateral<8> >()), m_element_rank);
174 register_cell_topology(fem::CellTopology(shards::getCellTopologyData< shards::ShellQuadrilateral<9> >()), m_element_rank);
180 ThrowRequireMsg(
is_FEM_initialized(),
"FEMMetaData::register_cell_topology: FEM_initialize() must be called before this function");
182 CellTopologyPartEntityRankMap::const_iterator it = m_cellTopologyPartEntityRankMap.find(cell_topology);
184 const bool duplicate = it != m_cellTopologyPartEntityRankMap.end();
185 const EntityRank existing_rank = duplicate ? (*it).second.second : 0;
187 ThrowInvalidArgMsgIf(m_spatial_dimension <
entity_rank,
189 "exceeds maximum spatial_dimension = " << m_spatial_dimension );
191 ThrowErrorMsgIf(duplicate && existing_rank !=
entity_rank,
192 "For args: cell_topolgy " << cell_topology.getName() <<
" and entity_rank " <<
entity_rank <<
", " <<
193 "previously declared rank = " << existing_rank );
196 std::string part_name = std::string(
"FEM_ROOT_CELL_TOPOLOGY_PART_") + std::string(cell_topology.getName());
198 ThrowErrorMsgIf(
get_part(part_name) != 0,
"Cannot register topology with same name as existing part '" << cell_topology.getName() <<
"'" );
201 m_cellTopologyPartEntityRankMap[cell_topology] = CellTopologyPartEntityRankMap::mapped_type(&part,
entity_rank);
203 assign_cell_topology(m_partCellTopologyVector, part.mesh_meta_data_ordinal(), cell_topology);
211 const std::string & topology_name)
const 213 std::string part_name = convert_to_internal_name(std::string(
"FEM_ROOT_CELL_TOPOLOGY_PART_") + topology_name);
219 return fem::CellTopology();
225 ThrowRequireMsg(
is_FEM_initialized(),
"FEMMetaData::get_cell_topology_root_part: FEM_initialize() must be called before this function");
226 CellTopologyPartEntityRankMap::const_iterator it = m_cellTopologyPartEntityRankMap.find(cell_topology);
227 ThrowErrorMsgIf(it == m_cellTopologyPartEntityRankMap.end(),
228 "Cell topology " << cell_topology.getName() <<
229 " has not been registered");
231 return *(*it).second.first;
241 ThrowRequireMsg(
is_FEM_initialized(),
"FEMMetaData::get_cell_topology: FEM_initialize() must be called before this function");
242 fem::CellTopology cell_topology;
245 if (part_ordinal < m_partCellTopologyVector.size())
247 cell_topology = m_partCellTopologyVector[part_ordinal];
250 return cell_topology;
254 void FEMMetaData::check_topo_db()
256 std::cout <<
"FEMMetaData::check_topo_db... m_partCellTopologyVector.size() = " << m_partCellTopologyVector.size() << std::endl;
258 fem::CellTopology cell_topology;
260 for (
unsigned i = 0; i < m_partCellTopologyVector.size(); i++)
262 cell_topology = m_partCellTopologyVector[i];
263 if (!cell_topology.getCellTopologyData())
265 std::cout <<
"bad topology in FEMMetaData::check_topo_db" << std::endl;
267 ThrowRequireMsg(cell_topology.getCellTopologyData(),
"bad topology in FEMMetaData::check_topo_db");
270 std::cout <<
"FEMMetaData::check_topo_db...done" << std::endl;
279 if (is_cell_topology_root_part(part)) {
283 for (PartVector::const_iterator it=subsets.begin() ; it != subsets.end() ; ++it) {
284 if (is_cell_topology_root_part( **it )) {
291 void find_cell_topologies_in_part_and_subsets_of_same_rank(
const Part & part, EntityRank rank, std::set<fem::CellTopology> & topologies_found)
294 fem::CellTopology top = fem_meta.get_cell_topology(part);
295 if ((top.isValid() && (part.primary_entity_rank() == rank))) {
296 topologies_found.insert(top);
299 for (PartVector::const_iterator it=subsets.begin() ; it != subsets.end() ; ++it) {
300 top = fem_meta.get_cell_topology(**it);
301 if (top.isValid() && ( (**it).primary_entity_rank() == rank) ) {
302 topologies_found.insert(top);
311 ThrowRequireMsg(
is_FEM_initialized(),
"FEMMetaData::declare_part_subset: FEM_initialize() must be called before this function");
314 const bool no_superset_topology = !superset_top.isValid();
315 if ( no_superset_topology ) {
320 const bool subset_has_root_part = root_part_in_subset(subset);
321 ThrowErrorMsgIf( subset_has_root_part,
"FEMMetaData::declare_part_subset: Error, root cell topology part found in subset or below." );
323 std::set<fem::CellTopology> cell_topologies;
324 find_cell_topologies_in_part_and_subsets_of_same_rank(subset,superset.
primary_entity_rank(),cell_topologies);
326 ThrowErrorMsgIf( cell_topologies.size() > 1,
327 "FEMMetaData::declare_part_subset: Error, multiple cell topologies of rank " 329 <<
" defined below subset" 331 const bool non_matching_cell_topology = ((cell_topologies.size() == 1) && (*cell_topologies.begin() != superset_top));
332 ThrowErrorMsgIf( non_matching_cell_topology,
333 "FEMMetaData::declare_part_subset: Error, superset topology = " 334 << superset_top.getName() <<
" does not match the topology = " 335 << cell_topologies.begin()->getName()
336 <<
" coming from the subset part" 344 for (PartVector::const_iterator it=subset_parts.begin() ; it != subset_parts.end() ; ++it) {
345 const Part & it_part = **it;
354 const fem::CellTopology cell_topology)
const 356 CellTopologyPartEntityRankMap::const_iterator it = m_cellTopologyPartEntityRankMap.find(cell_topology);
357 if (it == m_cellTopologyPartEntityRankMap.end())
360 return (*it).second.second;
367 bool is_cell_topology_root_part(
const Part & part) {
372 return (root_part == part);
382 void set_cell_topology(
384 fem::CellTopology cell_topology)
388 ThrowRequireMsg(fem_meta.is_FEM_initialized(),
"set_cell_topology: FEM_initialize() must be called before this function");
390 Part &root_part = fem_meta.get_cell_topology_root_part(cell_topology);
391 fem_meta.declare_part_subset(root_part, part);
394 std::vector<std::string>
395 entity_rank_names(
size_t spatial_dimension )
397 ThrowInvalidArgMsgIf( spatial_dimension < 1 || 3 < spatial_dimension,
398 "Invalid spatial dimension = " << spatial_dimension );
400 std::vector< std::string > names ;
402 names.reserve( spatial_dimension + 1 );
404 names.push_back( std::string(
"NODE" ) );
406 if ( 1 < spatial_dimension ) { names.push_back( std::string(
"EDGE") ); }
407 if ( 2 < spatial_dimension ) { names.push_back( std::string(
"FACE") ); }
409 names.push_back( std::string(
"ELEMENT") );
417 const Bucket & bucket)
419 const BulkData & bulk_data = BulkData::get(bucket);
420 const MetaData & meta_data = MetaData::get(bulk_data);
421 const PartVector & all_parts = meta_data.get_parts();
425 CellTopology cell_topology;
427 const std::pair< const unsigned *, const unsigned * > supersets = bucket.superset_part_ordinals();
429 if (supersets.first != supersets.second) {
430 const Part *first_found_part = 0;
432 for (
const unsigned * it = supersets.first ; it != supersets.second ; ++it ) {
434 const Part & part = * all_parts[*it] ;
436 if ( part.primary_entity_rank() == bucket.entity_rank() ) {
438 CellTopology top = fem.get_cell_topology( part );
440 if ( ! cell_topology.getCellTopologyData() ) {
441 cell_topology = top ;
443 if (!first_found_part)
444 first_found_part = ∂
447 ThrowErrorMsgIf( top.getCellTopologyData() && top != cell_topology,
448 "Cell topology is ambiguously defined. It is defined as " << cell_topology.getName() <<
449 " on part " << first_found_part->name() <<
" and as " << top.getName() <<
" on its superset part " << part.name() );
455 return cell_topology ;
unsigned primary_entity_rank() const
The primary entity type for this part.
An application-defined subset of a problem domain.
unsigned mesh_meta_data_ordinal() const
Internally generated ordinal of this part that is unique within the owning meta data manager...
const PartVector & subsets() const
Parts that are subsets of this part.
std::vector< Part *> PartVector
Collections of parts are frequently maintained as a vector of Part pointers.
EntityRank entity_rank(const EntityKey &key)
Given an entity key, return an entity type (rank).