Skip to content

Commit

Permalink
rename variable to be consistent with rest of class
Browse files Browse the repository at this point in the history
  • Loading branch information
white238 committed Nov 4, 2024
1 parent 784bc0c commit 35f54e6
Show file tree
Hide file tree
Showing 2 changed files with 23 additions and 13 deletions.
22 changes: 11 additions & 11 deletions src/serac/physics/state/state_manager.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -20,18 +20,18 @@ std::string StateManag
std::unordered_map<std::string, mfem::ParGridFunction*> StateManager::named_states_;
std::unordered_map<std::string, mfem::ParGridFunction*> StateManager::named_duals_;

double StateManager::newDataCollection(const std::string& name, const std::optional<int> cycle_to_load)
double StateManager::newDataCollection(const std::string& mesh_tag, const std::optional<int> cycle_to_load)
{
SLIC_ERROR_ROOT_IF(!ds_, "Cannot construct a DataCollection without a DataStore");
std::string coll_name = name + "_datacoll";
std::string coll_name = getCollectionName(mesh_tag);

auto global_grp = ds_->getRoot()->createGroup(coll_name + "_global");
auto bp_index_grp = global_grp->createGroup("blueprint_index/" + coll_name);
auto domain_grp = ds_->getRoot()->createGroup(coll_name);

// Needs to be configured to own the mesh data so all mesh data is saved to datastore/output file
constexpr bool owns_mesh_data = true;
auto [iter, _] = datacolls_.emplace(std::piecewise_construct, std::forward_as_tuple(name),
auto [iter, _] = datacolls_.emplace(std::piecewise_construct, std::forward_as_tuple(mesh_tag),
std::forward_as_tuple(coll_name, bp_index_grp, domain_grp, owns_mesh_data));
auto& datacoll = iter->second;
datacoll.SetComm(MPI_COMM_WORLD);
Expand All @@ -56,7 +56,7 @@ double StateManager::newDataCollection(const std::string& name, const std::optio
// indicates that the mesh is periodic and the new nodal grid function must also
// be discontinuous.
bool is_discontinuous = false;
auto nodes = mesh(name).GetNodes();
auto nodes = mesh(mesh_tag).GetNodes();
if (nodes) {
is_discontinuous = nodes->FESpace()->FEColl()->GetContType() == mfem::FiniteElementCollection::DISCONTINUOUS;
SLIC_WARNING_ROOT_IF(
Expand All @@ -71,17 +71,17 @@ double StateManager::newDataCollection(const std::string& name, const std::optio
// 2. Uses the existing continuity of the mesh finite element space (periodic meshes are discontinuous)
// 3. Uses the spatial dimension as the mesh dimension (i.e. it is not a lower dimension manifold)
// 4. Uses the ordering set by serac::ordering
mesh(name).SetCurvature(1, is_discontinuous, -1, serac::ordering);
mesh(mesh_tag).SetCurvature(1, is_discontinuous, -1, serac::ordering);

// Sidre will destruct the nodal grid function instead of the mesh
mesh(name).SetNodesOwner(false);
mesh(mesh_tag).SetNodesOwner(false);

// Generate the face neighbor information in the mesh. This is needed by the face restriction
// operators used by Functional
mesh(name).ExchangeFaceNbrData();
mesh(mesh_tag).ExchangeFaceNbrData();

// Construct and store the shape displacement fields and sensitivities associated with this mesh
constructShapeFields(name);
constructShapeFields(mesh_tag);

} else {
datacoll.SetCycle(0); // Iteration counter
Expand All @@ -95,9 +95,9 @@ void StateManager::loadCheckpointedStates(int cycle_to_load, std::vector<FiniteE
{
SERAC_MARK_FUNCTION;
mfem::ParMesh* meshPtr = &(*states_to_load.begin())->mesh();
std::string mesh_name = collectionID(meshPtr);
std::string mesh_tag = collectionID(meshPtr);

std::string coll_name = mesh_name + "_datacoll";
std::string coll_name = getCollectionName(mesh_tag);

axom::sidre::MFEMSidreDataCollection previous_datacoll(coll_name);

Expand All @@ -107,7 +107,7 @@ void StateManager::loadCheckpointedStates(int cycle_to_load, std::vector<FiniteE

for (auto state : states_to_load) {
meshPtr = &state->mesh();
SLIC_ERROR_ROOT_IF(collectionID(meshPtr) != mesh_name,
SLIC_ERROR_ROOT_IF(collectionID(meshPtr) != mesh_tag,
"Loading FiniteElementStates from two different meshes at one time is not allowed.");
mfem::ParGridFunction* datacoll_owned_grid_function = previous_datacoll.GetParField(state->name());

Expand Down
14 changes: 12 additions & 2 deletions src/serac/physics/state/state_manager.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -364,11 +364,21 @@ class StateManager {
private:
/**
* @brief Creates a new datacollection based on a registered mesh
* @param[in] name The name of the new datacollection
* @param[in] mesh_tag The mesh name used to name the new datacollection
* @param[in] cycle_to_load What cycle to load the DataCollection from, if applicable
* @return The time from specified restart cycle. Otherwise zero.
*/
static double newDataCollection(const std::string& name, const std::optional<int> cycle_to_load = {});
static double newDataCollection(const std::string& mesh_tag, const std::optional<int> cycle_to_load = {});

/**
* @brief Returns the name of the data collection's name for a given mesh_tag
* @param[in] mesh_tag The mesh name used to name the new datacollection
* @return The name of the data collection for the given mesh_tag
*/
static std::string getCollectionName(const std::string& mesh_tag)
{
return mesh_tag + "_datacoll";
}

/**
* @brief Construct the shape displacement field for the requested mesh
Expand Down

0 comments on commit 35f54e6

Please sign in to comment.