Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

OpenMCCellAverageProblem updates #907

Open
wants to merge 4 commits into
base: devel
Choose a base branch
from
Open

OpenMCCellAverageProblem updates #907

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
7715723
same changes to resolve conflicts again
Oct 20, 2024
1777287
clang-format
Oct 20, 2024
84329ed
updating tally system meshes
Oct 21, 2024
dd2b9df
clang-format
Oct 21, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
7 changes: 7 additions & 0 deletions include/base/OpenMCCellAverageProblem.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -21,6 +21,7 @@
#include "OpenMCProblemBase.h"
#include "SymmetryPointGenerator.h"
#include "OpenMCVolumeCalculation.h"
#include "MooseMesh.h"

/// Tally/filter includes.
#include "TallyBase.h"
Expand Down Expand Up @@ -870,6 +871,9 @@ class OpenMCCellAverageProblem : public OpenMCProblemBase
*/
const bool _needs_global_tally;

/// Whether OpenMCCellAverageProblem should use the displaced mesh
const bool & _use_displaced;

/**
* A map of the filter objects created by the [Problem/Filters] block. The key for each filter is
* it's corresponding MOOSE name to allow tallies to look up filters.
Expand Down Expand Up @@ -1030,6 +1034,9 @@ class OpenMCCellAverageProblem : public OpenMCProblemBase
/// Number of particles simulated in the first iteration
unsigned int _n_particles_1;

// Get a modifyable reference to the Moose mesh
virtual const MooseMesh & getMooseMesh() const;

/// Mapping from temperature variable name to the subdomains on which to read it from
std::map<std::string, std::vector<SubdomainName>> _temp_vars_to_blocks;

Expand Down
10 changes: 10 additions & 0 deletions include/tallies/CellTally.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -20,9 +20,13 @@

#include "TallyBase.h"
#include "OpenMCCellAverageProblem.h"
#include "MooseMesh.h"

#include "openmc/tallies/filter_cell.h"

class OpenMCCellAverageProblem;
class MooseMesh;

class CellTally : public TallyBase
{
public:
Expand Down Expand Up @@ -97,4 +101,10 @@ class CellTally : public TallyBase

/// Absolute tolerance for checking equal tally mapped volumes
const Real & _equal_tally_volume_abs_tol;

/// Whether the skinned mesh should be generated from a displaced mesh
const bool & _use_displaced;

/// Moose mesh
MooseMesh & getMooseMesh();
};
12 changes: 11 additions & 1 deletion include/tallies/MeshTally.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -20,9 +20,13 @@

#include "TallyBase.h"
#include "OpenMCCellAverageProblem.h"
#include "MooseMesh.h"

#include "openmc/tallies/filter_mesh.h"

class OpenMCCellAverageProblem;
class MooseMesh;

class MeshTally : public TallyBase
{
public:
Expand Down Expand Up @@ -63,7 +67,7 @@ class MeshTally : public TallyBase
* each element. This function performs as many checks as possible to ensure that the meshes
* are indeed identical.
*/
void checkMeshTemplateAndTranslations() const;
void checkMeshTemplateAndTranslations();

/**
* Mesh template file to use for creating mesh tallies in OpenMC; currently, this mesh
Expand All @@ -90,4 +94,10 @@ class MeshTally : public TallyBase

/// OpenMC unstructured mesh instance for use with mesh tallies
const openmc::LibMesh * _mesh_template;

/// Whether the skinned mesh should be generated from a displaced mesh
const bool & _use_displaced;

/// Moose mesh
MooseMesh & getMooseMesh();
};
7 changes: 7 additions & 0 deletions include/tallies/TallyBase.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,7 @@

#include "MooseObject.h"
#include "CardinalEnums.h"
#include "MooseMesh.h"

#include "openmc/tallies/tally.h"
#include "xtensor/xview.hpp"
Expand Down Expand Up @@ -301,4 +302,10 @@ class TallyBase : public MooseObject

/// Tolerance for setting zero tally
static constexpr Real ZERO_TALLY_THRESHOLD = 1e-12;

/// Whether the skinned mesh should be generated from a displaced mesh
const bool & _use_displaced;

/// Moose mesh
MooseMesh & getMooseMesh();
};
3 changes: 3 additions & 0 deletions include/userobjects/OpenMCVolumeCalculation.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -79,6 +79,9 @@ class OpenMCVolumeCalculation : public GeneralUserObject
/// Upper right of the box within which to compute OpenMC volumes
Point _upper_right;

/// Get a modifyable reference to the Moose mesh
MooseMesh & getMooseMesh();

/// Volume calculation object
std::unique_ptr<openmc::VolumeCalculation> _volume_calc;

Expand Down
66 changes: 43 additions & 23 deletions src/base/OpenMCCellAverageProblem.C
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,7 @@

#include "OpenMCCellAverageProblem.h"
#include "DelimitedFileReader.h"
#include "DisplacedProblem.h"
#include "TallyBase.h"
#include "CellTally.h"
#include "AddTallyAction.h"
Expand Down Expand Up @@ -164,7 +165,9 @@ OpenMCCellAverageProblem::validParams()
params.addParam<int>("first_iteration_particles",
"Number of particles to use for first iteration "
"when using Dufek-Gudowski relaxation");

params.addParam<bool>("use_displaced_mesh",
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@nuclearkevin does this need to be modified at all to be consistent with what you did with fixed_mesh being auto-detected?

@meltawila can you please also remind me, was there a use case where we wanted to disable the displaced mesh updates in OpenMC, even if the mesh is actually moving? Or was this user option only existing because we didn't know how to detect it automatically?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is the logic used in MOOSE solid mechanics, to be able to set use_displaced_mesh = false for certain post processors/user objects only for example

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@aprilnovak At the moment this doesn't need to be modified since the changes made to auto-detect fixed_mesh are coming in the AMR PR (#964). If this gets merged first, I can handle the necessary changes in #964. If the AMR PR gets merged ahead of this one, a single line of code will need to be added to this PR to get auto-detection working:
https://github.com/neams-th-coe/cardinal/pull/964/files#diff-7321f62ca7502f529d9917464fed4643c5e57c97c27d295f78db837d1a35a366R219

Copy link
Contributor

@nuclearkevin nuclearkevin Oct 20, 2024 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

For reference @meltawila this is what's necessary to automatically detect if the mesh is moving and setup the OpenMCCellAverageProblem accordingly. This include is required:

#include "CreateDisplacedProblemAction.h"

and this needs to get added to the constructor of OpenMCCellAverageProblem:

// Check to see if a displaced problem is being initialized.
const auto & dis_actions = getMooseApp().actionWarehouse().getActions<CreateDisplacedProblemAction>();
for (const auto & act : dis_actions)
{
  auto has_displaced = act->isParamValid("displacements") && act->getParam<bool>("use_displaced_mesh");
  _need_to_reinit_coupling |= (has_displaced && _use_displaced_mesh);
}

I already have this in a separate PR, so if yours gets merged first you won't need to add this.

false,
"Whether OpenMCCellAverageProblem should use the displaced mesh ");
params.addParam<UserObjectName>(
"symmetry_mapper",
"User object (of type SymmetryPointGenerator) "
Expand Down Expand Up @@ -207,6 +210,7 @@ OpenMCCellAverageProblem::OpenMCCellAverageProblem(const InputParameters & param
_specified_temperature_feedback(params.isParamSetByUser("temperature_blocks")),
_needs_to_map_cells(_specified_density_feedback || _specified_temperature_feedback),
_needs_global_tally(_check_tally_sum || _normalize_by_global),
_use_displaced(getParam<bool>("use_displaced_mesh")),
_volume_calc(nullptr),
_symmetry(nullptr)
{
Expand Down Expand Up @@ -368,6 +372,16 @@ OpenMCCellAverageProblem::OpenMCCellAverageProblem(const InputParameters & param
}
}

const MooseMesh &
OpenMCCellAverageProblem::getMooseMesh() const
{
if (_use_displaced && !_displaced_problem && !_first_transfer)
mooseWarning("Displaced mesh was requested but the displaced problem does not exist. "
"Regular mesh will be returned");
return ((_use_displaced && _displaced_problem && !_first_transfer) ? _displaced_problem->mesh()
: mesh());
}

void
OpenMCCellAverageProblem::readBlockVariables(
const std::string & param,
Expand Down Expand Up @@ -545,9 +559,9 @@ OpenMCCellAverageProblem::setupProblem()
{
// establish the local -> global element mapping for convenience
_local_to_global_elem.clear();
for (unsigned int e = 0; e < _mesh.nElem(); ++e)
for (unsigned int e = 0; e < getMooseMesh().nElem(); ++e)
{
const auto * elem = _mesh.queryElemPtr(e);
const auto * elem = getMooseMesh().queryElemPtr(e);
if (!isLocalElem(elem))
continue;

Expand Down Expand Up @@ -633,7 +647,7 @@ OpenMCCellAverageProblem::readBlockParameters(const std::string name,
auto names = getParam<std::vector<SubdomainName>>(name);
checkEmptyVector(names, "'" + name + "'");

auto b_ids = _mesh.getSubdomainIDs(names);
auto b_ids = getMooseMesh().getSubdomainIDs(names);
std::copy(b_ids.begin(), b_ids.end(), std::inserter(blocks, blocks.end()));
checkBlocksInMesh(name, b_ids, names);
}
Expand All @@ -644,7 +658,7 @@ OpenMCCellAverageProblem::checkBlocksInMesh(const std::string name,
const std::vector<SubdomainID> & ids,
const std::vector<SubdomainName> & names) const
{
const auto & subdomains = _mesh.meshSubdomains();
const auto & subdomains = getMooseMesh().meshSubdomains();
for (std::size_t b = 0; b < names.size(); ++b)
if (subdomains.find(ids[b]) == subdomains.end())
mooseError("Block '" + names[b] + "' specified in '" + name + "' " + "not found in mesh!");
Expand Down Expand Up @@ -672,7 +686,7 @@ OpenMCCellAverageProblem::read2DBlockParameters(const std::string name,
for (const auto & i : slice)
flattened_names.push_back(i);

flattened_ids = _mesh.getSubdomainIDs(flattened_names);
flattened_ids = getMooseMesh().getSubdomainIDs(flattened_names);
checkBlocksInMesh(name, flattened_ids, flattened_names);

// should not be any duplicate blocks
Expand Down Expand Up @@ -742,8 +756,8 @@ OpenMCCellAverageProblem::storeElementPhase()
for (const auto & s : excl_density_blocks)
_n_moose_density_elems += numElemsInSubdomain(s);

_n_moose_none_elems =
_mesh.nElem() - _n_moose_temp_density_elems - _n_moose_temp_elems - _n_moose_density_elems;
_n_moose_none_elems = getMooseMesh().nElem() - _n_moose_temp_density_elems - _n_moose_temp_elems -
_n_moose_density_elems;
}

void
Expand All @@ -757,7 +771,7 @@ OpenMCCellAverageProblem::computeCellMappedVolumes()
for (const auto & e : c.second)
{
// we are looping over local elements, so no need to check for nullptr
const auto * elem = _mesh.queryElemPtr(globalElemID(e));
const auto * elem = getMooseMesh().queryElemPtr(globalElemID(e));
vol += elem->volume();
}

Expand Down Expand Up @@ -837,7 +851,7 @@ OpenMCCellAverageProblem::getCellMappedPhase()

// we are looping over local elements, so no need to check for nullptr
for (const auto & e : c.second)
f[elemFeedback(_mesh.queryElemPtr(globalElemID(e)))]++;
f[elemFeedback(getMooseMesh().queryElemPtr(globalElemID(e)))]++;

cells_n_temp.push_back(f[coupling::temperature]);
cells_n_temp_rho.push_back(f[coupling::density_and_temperature]);
Expand Down Expand Up @@ -998,7 +1012,7 @@ OpenMCCellAverageProblem::printAuxVariableIO()
VariadicTable<std::string, std::string, std::string> aux(
{"Subdomain", "Temperature", "Density"});

for (const auto & s : _mesh.meshSubdomains())
for (const auto & s : getMooseMesh().meshSubdomains())
{
std::string temp = _subdomain_to_temp_vars.count(s) ? _subdomain_to_temp_vars[s].second : "";
std::string rho =
Expand Down Expand Up @@ -1058,7 +1072,7 @@ OpenMCCellAverageProblem::getCellMappedSubdomains()
for (const auto & e : c.second)
{
// we are looping over local elements, so no need to check for nullptr
const auto * elem = _mesh.queryElemPtr(globalElemID(e));
const auto * elem = getMooseMesh().queryElemPtr(globalElemID(e));
elem_ids.push_back(elem->subdomain_id());
}
}
Expand Down Expand Up @@ -1325,9 +1339,10 @@ OpenMCCellAverageProblem::initializeElementToCellMapping()
mooseError("Did not find any overlap between MOOSE elements and OpenMC cells for "
"the specified blocks!");

_console << "\nMapping between " + Moose::stringify(_mesh.nElem()) + " MOOSE elements and " +
Moose::stringify(_n_openmc_cells) + " OpenMC cells (on " +
Moose::stringify(openmc::model::n_coord_levels) + " coordinate levels):"
_console << "\nMapping between " + Moose::stringify(getMooseMesh().nElem()) +
" MOOSE elements and " + Moose::stringify(_n_openmc_cells) +
" OpenMC cells (on " + Moose::stringify(openmc::model::n_coord_levels) +
" coordinate levels):"
<< std::endl;

VariadicTable<std::string, int, int, int, int> vt(
Expand Down Expand Up @@ -1633,13 +1648,14 @@ OpenMCCellAverageProblem::mapElemsToCells()
// reset data structures
_elem_to_cell.clear();
_cell_to_elem.clear();
_subdomain_to_material.clear();
_flattened_ids.clear();
_flattened_instances.clear();

int local_elem = -1;
for (unsigned int e = 0; e < _mesh.nElem(); ++e)
for (unsigned int e = 0; e < getMooseMesh().nElem(); ++e)
{
const auto * elem = _mesh.queryElemPtr(e);
const auto * elem = getMooseMesh().queryElemPtr(e);

if (!isLocalElem(elem))
continue;
Expand Down Expand Up @@ -1764,8 +1780,8 @@ OpenMCCellAverageProblem::mapElemsToCells()
gatherCellVector(elems, n_elems, _cell_to_elem);

// fill out the elem_to_cell structure
_elem_to_cell.resize(_mesh.nElem());
for (unsigned int e = 0; e < _mesh.nElem(); ++e)
_elem_to_cell.resize(getMooseMesh().nElem());
for (unsigned int e = 0; e < getMooseMesh().nElem(); ++e)
_elem_to_cell[e] = {UNMAPPED, UNMAPPED};

for (const auto & c : _cell_to_elem)
Expand All @@ -1784,7 +1800,7 @@ OpenMCCellAverageProblem::getPointInCell()
for (const auto & c : _local_cell_to_elem)
{
// we are only dealing with local elements here, no need to check for nullptr
const Elem * elem = _mesh.queryElemPtr(globalElemID(c.second[0]));
const Elem * elem = getMooseMesh().queryElemPtr(globalElemID(c.second[0]));
const Point & p = elem->vertex_average();

x.push_back(p(0));
Expand Down Expand Up @@ -1974,7 +1990,7 @@ OpenMCCellAverageProblem::addExternalVariables()
{
auto number = addExternalVariable(v.first, &v.second);

auto ids = _mesh.getSubdomainIDs(v.second);
auto ids = getMooseMesh().getSubdomainIDs(v.second);
for (const auto & s : ids)
_subdomain_to_density_vars[s] = {number, v.first};
}
Expand All @@ -1984,7 +2000,7 @@ OpenMCCellAverageProblem::addExternalVariables()
{
auto number = addExternalVariable(v.first, &v.second);

auto ids = _mesh.getSubdomainIDs(v.second);
auto ids = getMooseMesh().getSubdomainIDs(v.second);
for (const auto & s : ids)
_subdomain_to_temp_vars[s] = {number, v.first};
}
Expand Down Expand Up @@ -2049,7 +2065,7 @@ OpenMCCellAverageProblem::computeVolumeWeightedCellInput(
for (const auto & e : c.second)
{
// we are only accessing local elements here, so no need to check for nullptr
const auto * elem = _mesh.queryElemPtr(globalElemID(e));
const auto * elem = getMooseMesh().queryElemPtr(globalElemID(e));
auto v = var_num.at(elem->subdomain_id()).first;
auto dof_idx = elem->dof_number(sys_number, v, 0);
product += (*_serialized_solution)(dof_idx) * elem->volume();
Expand Down Expand Up @@ -2322,6 +2338,10 @@ OpenMCCellAverageProblem::syncSolutions(ExternalProblem::Direction direction)
if (_volume_calc)
_volume_calc->resetVolumeCalculation();

if (_use_displaced)
meltawila marked this conversation as resolved.
Show resolved Hide resolved
{
_displaced_problem->updateMesh();
}
resetTallies();
setupProblem();
}
Expand Down
Loading