43#include "PanzerAdaptersSTK_config.hpp"
47#include "Teuchos_TimeMonitor.hpp"
48#include "Teuchos_StandardParameterEntryValidators.hpp"
49#include <stk_mesh/base/FieldBase.hpp>
59 stk::ParallelMachine mpi_comm,
60 const bool print_debug)
61 : createEdgeBlocks_(false),
62 print_debug_(print_debug)
64 panzer_stk::STK_ExodusReaderFactory factory;
65 RCP<Teuchos::ParameterList> pl = rcp(
new Teuchos::ParameterList);
66 pl->set(
"File Name",quad8MeshFileName);
67 factory.setParameterList(pl);
74 const bool print_debug)
75 : quad8Mesh_(quad8Mesh),
76 createEdgeBlocks_(false),
77 print_debug_(print_debug)
85 PANZER_FUNC_TIME_MONITOR(
"panzer::Quad8ToQuad4MeshFactory::buildMesh()");
90 MPI_Comm_compare(parallelMach,
quad8Mesh_->getBulkData()->parallel(), &
result);
91 TEUCHOS_ASSERT(
result != MPI_UNEQUAL);
98 mesh->initialize(parallelMach);
108 PANZER_FUNC_TIME_MONITOR(
"panzer::Quad8ToQuad4MeshFactory::buildUncomittedMesh()");
112 machRank_ = stk::parallel_machine_rank(parallelMach);
113 machSize_ = stk::parallel_machine_size(parallelMach);
126 PANZER_FUNC_TIME_MONITOR(
"panzer::Quad8ToQuad4MeshFactory::completeMeshConstruction()");
135#ifndef ENABLE_UNIFORM
144#ifndef ENABLE_UNIFORM
167 setMyParamList(paramList);
180 static RCP<Teuchos::ParameterList> defaultParams;
183 if(defaultParams == Teuchos::null) {
184 defaultParams = rcp(
new Teuchos::ParameterList);
186 Teuchos::setStringToIntegralParameter<int>(
187 "Offset mesh GIDs above 32-bit int limit",
189 "If 64-bit GIDs are supported, the mesh element and node global indices will start at a value greater than 32-bit limit.",
190 Teuchos::tuple<std::string>(
"OFF",
"ON"),
191 defaultParams.get());
194 defaultParams->set<
bool>(
"Create Edge Blocks",
false,
"Create edge blocks in the mesh");
196 Teuchos::ParameterList & bcs = defaultParams->sublist(
"Periodic BCs");
197 bcs.set<
int>(
"Count",0);
200 return defaultParams;
206 std::cout <<
"\n\n**** DEBUG: begin printing source Quad8 exodus file metadata ****\n";
207 quad8Mesh_->getMetaData()->dump_all_meta_info(std::cout);
208 std::cout <<
"\n\n**** DEBUG: end printing source Quad8 exodus file metadata ****\n";
212 using QuadTopo = shards::Quadrilateral<4>;
213 const CellTopologyData * ctd = shards::getCellTopologyData<QuadTopo>();
214 const CellTopologyData * side_ctd = shards::CellTopology(ctd).getBaseCellTopologyData(1,0);
217 std::vector<std::string> element_block_names;
218 quad8Mesh_->getElementBlockNames(element_block_names);
220 for (
const auto& n : element_block_names)
226 std::vector<std::string> sideset_names;
228 for (
const auto& n : sideset_names)
234 std::vector<std::string> nodeset_names;
236 for (
const auto& n : nodeset_names)
241 const CellTopologyData * edge_ctd = shards::CellTopology(ctd).getBaseCellTopologyData(1,0);
242 std::vector<std::string> element_block_names;
243 quad8Mesh_->getElementBlockNames(element_block_names);
244 for (
const auto& block_name : element_block_names)
251 for (
const auto& f : fields) {
253 std::cout <<
"Field=" << f->name() <<
", rank=" << f->entity_rank() << std::endl;
257 if (f->name() !=
"coordinates") {
258 for (
const auto& n : element_block_names)
267 for (
const auto& f : fields) {
269 std::cout <<
"Add Cell Field=" << f->name() <<
", rank=" << f->entity_rank() << std::endl;
271 for (
const auto& n : element_block_names)
283 std::cout <<
"\n\n**** DEBUG: begin printing source Quad4 exodus file metadata ****\n";
285 std::cout <<
"\n\n**** DEBUG: end printing source Quad4 exodus file metadata ****\n";
297 std::vector<std::string> block_names;
298 quad8Mesh_->getElementBlockNames(block_names);
299 for (
const auto& block_name : block_names) {
302 std::vector<stk::mesh::Entity> elements;
303 quad8Mesh_->getMyElements(block_name,elements);
306 std::cout <<
"*************************************************" << std::endl;
307 std::cout <<
"block_name=" << block_name <<
", num_my_elements=" << elements.size() << std::endl;
308 std::cout <<
"*************************************************" << std::endl;
311 for (
const auto& element : elements) {
313 const auto element_gid =
quad8Mesh_->getBulkData()->identifier(element);
317 <<
", block=" << block_name
318 <<
", element entity_id=" << element
319 <<
", gid=" << element_gid << std::endl;
323 std::vector<stk::mesh::EntityId> nodes(4);
324 for (
int i=0; i < 4; ++i) {
331 TEUCHOS_ASSERT(node_entity.is_local_offset_valid());
332 const auto node_gid =
quad8Mesh_->getBulkData()->identifier(node_entity);
333 const double* node_coords =
quad8Mesh_->getNodeCoordinates(node_entity);
334 std::vector<double> coords_vec(2);
335 coords_vec[0] = node_coords[0];
336 coords_vec[1] = node_coords[1];
337 mesh.
addNode(node_gid,coords_vec);
340 std::cout <<
"elem gid=" <<
quad8Mesh_->getBulkData()->identifier(element)
341 <<
", node_gid=" << node_gid <<
", (" << coords_vec[0] <<
"," << coords_vec[1] <<
")\n";
347 auto element_block_part = mesh.
getMetaData()->get_part(block_name);
348 mesh.
addElement(element_descriptor,element_block_part);
359 std::vector<std::string> sideset_names;
361 for (
const auto& sideset_name : sideset_names) {
363 stk::mesh::Part* sideset_part = mesh.
getSideset(sideset_name);
365 std::vector<stk::mesh::Entity> q8_sides;
366 quad8Mesh_->getMySides(sideset_name,q8_sides);
369 for (
const auto q8_edge : q8_sides) {
374 stk::mesh::EntityId edge_gid =
quad8Mesh_->getBulkData()->identifier(q8_edge);
390 Teuchos::RCP<stk::mesh::BulkData> bulkData = mesh.
getBulkData();
391 Teuchos::RCP<stk::mesh::MetaData> metaData = mesh.
getMetaData();
395 stk::mesh::Selector owned_block = metaData->locally_owned_part();
397 std::vector<stk::mesh::Entity> edges;
398 bulkData->get_entities(mesh.
getEdgeRank(), owned_block, edges);
410 for (
const auto&
field : fields) {
413 std::cout <<
"Adding field values for \"" << *
field <<
"\" to the Quad4 mesh!\n";
415 auto field_name =
field->name();
419 if (
field->type_is<
double>() &&
420 field_name !=
"coordinates" &&
421 field_name !=
"PROC_ID" &&
422 field_name !=
"LOAD_BAL") {
425 std::vector<std::string> block_names;
426 quad8Mesh_->getElementBlockNames(block_names);
427 for (
const auto& block : block_names) {
429 auto* q4_field = q4_mesh.
getCellField(field_name,block);
430 TEUCHOS_ASSERT(q4_field !=
nullptr);
435 TEUCHOS_ASSERT(q8_field !=
nullptr);
439 std::vector<stk::mesh::Entity> q4_elements;
441 std::vector<stk::mesh::Entity> q8_elements;
443 TEUCHOS_ASSERT(q4_elements.size() == q8_elements.size());
445 for (
size_t i=0; i < q4_elements.size(); ++i) {
447 TEUCHOS_ASSERT(q4_mesh.
getBulkData()->identifier(q4_elements[i]) ==
448 quad8Mesh_->getBulkData()->identifier(q8_elements[i]));
451 double*
const q4_val =
static_cast<double*
>(stk::mesh::field_data(*q4_field,q4_elements[i]));
452 const double*
const q8_val =
static_cast<double*
>(stk::mesh::field_data(*q8_field,q8_elements[i]));
456 std::cout <<
"field=" << field_name <<
", block=" << block
457 <<
", q4e(" << q4_mesh.
getBulkData()->identifier(q4_elements[i]) <<
") = " << *q4_val
458 <<
", q8e(" <<
quad8Mesh_->getBulkData()->identifier(q8_elements[i]) <<
") = " << *q8_val
PHX::MDField< ScalarT, panzer::Cell, panzer::IP > result
A field that will be used to build up the result of the integral we're performing.
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.
virtual void completeMeshConstruction(STK_Interface &mesh, stk::ParallelMachine parallelMach) const
void addSideSets(STK_Interface &mesh) const
void addNodeSets(STK_Interface &mesh) const
Quad8ToQuad4MeshFactory(const std::string &quad8MeshFileName, stk::ParallelMachine mpi_comm=MPI_COMM_WORLD, const bool print_debug=false)
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
Derived from ParameterListAcceptor.
Teuchos::RCP< STK_Interface > buildMesh(stk::ParallelMachine parallelMach) const
Build the mesh object.
void buildElements(stk::ParallelMachine parallelMach, STK_Interface &mesh) const
void buildMetaData(stk::ParallelMachine parallelMach, STK_Interface &mesh) const
Teuchos::RCP< panzer_stk::STK_Interface > quad8Mesh_
void addEdgeBlocks(STK_Interface &mesh) const
virtual Teuchos::RCP< STK_Interface > buildUncommitedMesh(stk::ParallelMachine parallelMach) const
void copyCellFieldData(STK_Interface &mesh) const
std::string edgeBlockName_
void setParameterList(const Teuchos::RCP< Teuchos::ParameterList > ¶mList)
Derived from ParameterListAcceptor.
void buildLocalElementIDs()
void initialize(stk::ParallelMachine parallelMach, bool setupIO=true, const bool buildRefinementSupport=false)
static const std::string edgeBlockString
stk::mesh::Part * getEdgeBlock(const std::string &name) const
get the block part
bool isInitialized() const
Has initialize been called on this mesh object?
stk::mesh::EntityRank getNodeRank() const
void addEdgeBlock(const std::string &elemBlockName, const std::string &edgeBlockName, const stk::topology &topology)
void buildSubcells()
force the mesh to build subcells: edges and faces
stk::mesh::EntityRank getElementRank() const
void addElement(const Teuchos::RCP< ElementDescriptor > &ed, stk::mesh::Part *block)
void addNode(stk::mesh::EntityId gid, const std::vector< double > &coord)
void addCellField(const std::string &fieldName, const std::string &blockId)
void addNodeset(const std::string &name)
void addSideset(const std::string &name, const CellTopologyData *ctData)
Teuchos::RCP< stk::mesh::MetaData > getMetaData() const
stk::mesh::Field< double > * getCellField(const std::string &fieldName, const std::string &blockId) const
void addEntitiesToEdgeBlock(std::vector< stk::mesh::Entity > entities, stk::mesh::Part *edgeblock)
void addSolutionField(const std::string &fieldName, const std::string &blockId)
void addEntityToSideset(stk::mesh::Entity entity, stk::mesh::Part *sideset)
void getMyElements(std::vector< stk::mesh::Entity > &elements) const
void addElementBlock(const std::string &name, const CellTopologyData *ctData)
Teuchos::RCP< stk::mesh::BulkData > getBulkData() const
stk::mesh::EntityRank getEdgeRank() const
stk::mesh::Part * getSideset(const std::string &name) const
static void parsePeriodicBCList(const Teuchos::RCP< Teuchos::ParameterList > &pl, std::vector< Teuchos::RCP< const PeriodicBC_MatcherBase > > &periodicBC, bool &useBBoxSearch)
std::vector< Teuchos::RCP< const PeriodicBC_MatcherBase > > periodicBCVec_
Teuchos::RCP< ElementDescriptor > buildElementDescriptor(stk::mesh::EntityId elmtId, std::vector< stk::mesh::EntityId > &nodes)