44#include <Teuchos_TimeMonitor.hpp>
45#include <PanzerAdaptersSTK_config.hpp>
63 Teuchos::RCP<STK_Interface>
66 PANZER_FUNC_TIME_MONITOR(
"panzer::CustomMeshFactory::buildMesh()");
72 mesh->initialize(parallelMach);
80 Teuchos::RCP<STK_Interface>
83 PANZER_FUNC_TIME_MONITOR(
"panzer::CustomMeshFactory::buildUncomittedMesh()");
87 machRank_ = stk::parallel_machine_rank(parallelMach);
88 machSize_ = stk::parallel_machine_size(parallelMach);
101 stk::ParallelMachine parallelMach)
const
103 PANZER_FUNC_TIME_MONITOR(
"panzer::CustomMeshFactory::completeMeshConstruction()");
133 setMyParamList(paramList);
135 Dimension_ = paramList->get<
int>(
"Dimension");
137 NumBlocks_ = paramList->get<
int>(
"NumBlocks");
140 Nodes_ = paramList->get<
int*>(
"Nodes");
142 Coords_ = paramList->get<
double*>(
"Coords");
146 BlockIDs_ = paramList->get<
int*>(
"BlockIDs");
160 Teuchos::RCP<const Teuchos::ParameterList>
163 static RCP<Teuchos::ParameterList> defaultParams;
166 if(defaultParams == Teuchos::null) {
167 defaultParams = rcp(
new Teuchos::ParameterList);
169 defaultParams->set<
int>(
"Dimension",3);
171 defaultParams->set<
int>(
"NumBlocks",0);
173 defaultParams->set<
int>(
"NumNodesPerProc",0);
174 defaultParams->set<
int*>(
"Nodes",NULL);
176 defaultParams->set<
double*>(
"Coords",NULL);
178 defaultParams->set<
int>(
"NumElementsPerProc",0);
180 defaultParams->set<
int*>(
"BlockIDs",NULL);
181 defaultParams->set<
int*>(
"Element2Nodes",NULL);
183 defaultParams->set<
int>(
"OffsetToGlobalElementIDs", 0);
185 defaultParams->set<
double*>(
"ChargeDensity",NULL);
186 defaultParams->set<
double*>(
"ElectricPotential",NULL);
188 Teuchos::ParameterList &bcs = defaultParams->sublist(
"Periodic BCs");
189 bcs.set<
int>(
"Count",0);
192 return defaultParams;
199 RCP<Teuchos::ParameterList> validParams = rcp(
new Teuchos::ParameterList(*
getValidParameters()));
208 typedef shards::Hexahedron<8> HexTopo;
209 const CellTopologyData *ctd = shards::getCellTopologyData<HexTopo>();
210 const CellTopologyData *side_ctd = shards::CellTopology(ctd).getBaseCellTopologyData(2,0);
213 std::stringstream block_id;
214 block_id <<
"eblock-" << blk;
241 std::vector<double> coords(dim,0.0);
243 for (
int k=0;k<dim;++k)
249 std::vector<std::string> block_ids;
255 std::stringstream block_id;
262 std::vector<stk::mesh::EntityId> elt2nodes(8);
264 for (
int k=0;k<8;++k)
278 stk::mesh::BulkData& bulkData = *mesh.
getBulkData();
279 const stk::mesh::EntityRank sideRank = mesh.
getSideRank();
282 stk::mesh::Part *box[6];
290 stk::mesh::Part *wall = mesh.
getSideset(
"wall");
292 std::vector<stk::mesh::Entity> elements;
296 for (std::vector<stk::mesh::Entity>::const_iterator
297 itr=elements.begin();itr!=elements.end();++itr) {
298 stk::mesh::Entity element = (*itr);
299 const size_t numSides = bulkData.num_connectivity(element, sideRank);
300 stk::mesh::Entity
const* relations = bulkData.begin(element, sideRank);
303 for (std::size_t i=0;i<numSides;++i) {
304 stk::mesh::Entity side = relations[i];
305 const size_t numNeighbors = bulkData.num_elements(side);
306 stk::mesh::Entity
const* neighbors = bulkData.begin_elements(side);
308 if (numNeighbors == 1) {
312 else if (numNeighbors == 2) {
330 constexpr std::size_t dim_1 = 8;
334 std::stringstream block_id;
335 block_id <<
"eblock-" << blk;
338 std::vector<stk::mesh::Entity> elements;
342 const auto n_elements = elements.size();
345 std::vector<std::size_t> local_ids;
347 Kokkos::View<double**, Kokkos::HostSpace> charge_density_by_local_ids(
"charge_density_by_local_ids", n_elements, dim_1);
348 Kokkos::View<double**, Kokkos::HostSpace> electric_potential_by_local_ids(
"electric_potential_by_local_ids", n_elements, dim_1);
350 for (
const auto& elem : elements){
354 for (std::size_t k=0;k<dim_1;++k) {
355 const auto loc = q*dim_1 + k;
365 charge_density_by_local_ids, 1.0);
370 electric_potential_by_local_ids, 1.0);
double * ElectricPotential_
void buildMetaData(STK_Interface &mesh) const
virtual ~CustomMeshFactory()
Destructor.
void addSideSets(STK_Interface &mesh) const
int OffsetToGlobalElementIDs_
void fillSolutionFieldData(STK_Interface &mesh) const
void setParameterList(const Teuchos::RCP< Teuchos::ParameterList > ¶mList)
From ParameterListAcceptor.
CustomMeshFactory()
Constructor.
void initializeWithDefaults()
virtual void completeMeshConstruction(STK_Interface &mesh, stk::ParallelMachine parallelMach) const
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
From ParameterListAcceptor.
virtual Teuchos::RCP< STK_Interface > buildUncommitedMesh(stk::ParallelMachine parallelMach) const
void buildElements(STK_Interface &mesh) const
Teuchos::RCP< STK_Interface > buildMesh(stk::ParallelMachine parallelMach) const
Build the mesh object.
void buildLocalElementIDs()
void initialize(stk::ParallelMachine parallelMach, bool setupIO=true, const bool buildRefinementSupport=false)
stk::mesh::Part * getElementBlockPart(const std::string &name) const
get the block part
stk::mesh::EntityId elementGlobalId(std::size_t lid) const
bool isInitialized() const
Has initialize been called on this mesh object?
std::string containingBlockId(stk::mesh::Entity elmt) const
void buildSubcells()
force the mesh to build subcells: edges and faces
void addElement(const Teuchos::RCP< ElementDescriptor > &ed, stk::mesh::Part *block)
void addNode(stk::mesh::EntityId gid, const std::vector< double > &coord)
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::size_t elementLocalId(stk::mesh::Entity elmt) const
void getElementBlockNames(std::vector< std::string > &names) const
void addSideset(const std::string &name, const CellTopologyData *ctData)
unsigned getDimension() const
get the dimension
unsigned entityOwnerRank(stk::mesh::Entity entity) const
void addSolutionField(const std::string &fieldName, const std::string &blockId)
void addEntityToSideset(stk::mesh::Entity entity, stk::mesh::Part *sideset)
stk::mesh::EntityRank getSideRank() const
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::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_