Skip to content

Commit

Permalink
adding reporting ID for plane, refs idaholab#19217
Browse files Browse the repository at this point in the history
  • Loading branch information
yjung-anl committed Dec 2, 2021
1 parent ece47c7 commit 0a58ec5
Show file tree
Hide file tree
Showing 9 changed files with 279 additions and 0 deletions.
37 changes: 37 additions & 0 deletions modules/reactor/include/meshgenerators/PlaneIDMeshGenerator.h
@@ -0,0 +1,37 @@
//* This file is part of the MOOSE framework
//* https://www.mooseframework.org
//*
//* All rights reserved, see COPYRIGHT for full restrictions
//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
//*
//* Licensed under LGPL 2.1, please see LICENSE for details
//* https://www.gnu.org/licenses/lgpl-2.1.html

#pragma once

#include "MeshGenerator.h"

/**
* Assigns plane extra ID for existing 3D meshes
*/
class PlaneIDMeshGenerator : public MeshGenerator
{
public:
static InputParameters validParams();

PlaneIDMeshGenerator(const InputParameters & parameters);

std::unique_ptr<MeshBase> generate() override;

protected:
/// get plane ID for given plane coordiantes
int getPlaneID(const Point & p) const;
/// input mesh for adding reporting ID
std::unique_ptr<MeshBase> & _input;
/// coordinates of planes
std::vector<Real> _planes;
/// index of plane axis
const unsigned int _axis_index;
/// name of integer ID
const std::string _element_id_name;
};
110 changes: 110 additions & 0 deletions modules/reactor/src/meshgenerators/PlaneIDMeshGenerator.C
@@ -0,0 +1,110 @@
//* This file is part of the MOOSE framework
//* https://www.mooseframework.org
//*
//* All rights reserved, see COPYRIGHT for full restrictions
//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
//*
//* Licensed under LGPL 2.1, please see LICENSE for details
//* https://www.gnu.org/licenses/lgpl-2.1.html

#include "PlaneIDMeshGenerator.h"
#include "CastUniquePointer.h"

#include "libmesh/elem.h"

registerMooseObject("MooseApp", PlaneIDMeshGenerator);

InputParameters
PlaneIDMeshGenerator::validParams()
{
InputParameters params = MeshGenerator::validParams();

params.addRequiredParam<MeshGeneratorName>("input", "The mesh we want to modify");
params.addRequiredParam<std::vector<Real>>("plane_coordinates", "Coordinates of planes");
params.addParam<std::vector<unsigned int>>("num_ids_per_plane", "Number of unique ids per plane");
MooseEnum plane_axis("x y z", "z");
params.addParam<MooseEnum>("plane_axis", plane_axis, "Axis of plane");
params.addRequiredParam<std::string>("id_name", "Name of extra integer ID set");
params.addParam<Real>("tolerance", 1.0E-4, "Tolerance for plane coordinate check");
params.addClassDescription("Adds an extra element integer that identifies planes in a mesh.");
return params;
}

PlaneIDMeshGenerator::PlaneIDMeshGenerator(const InputParameters & parameters)
: MeshGenerator(parameters),
_input(getMesh("input")),
_axis_index(getParam<MooseEnum>("plane_axis")),
_element_id_name(getParam<std::string>("id_name"))
{
if (!isParamValid("num_ids_per_plane"))
_planes = getParam<std::vector<Real>>("plane_coordinates");
else
{
_planes.clear();
std::vector<Real> base_planes = getParam<std::vector<Real>>("plane_coordinates");
std::vector<unsigned int> sublayers = getParam<std::vector<unsigned int>>("num_ids_per_plane");
if (base_planes.size() != sublayers.size() + 1)
paramError("plane_coordinates",
"Sizes of 'plane_coordinates' and 'num_ids_per_plane' disagree");
_planes.push_back(base_planes[0]);
for (unsigned int i = 0; i < sublayers.size(); ++i)
{
Real layer_size = (base_planes[i + 1] - base_planes[i]) / (Real)sublayers[i];
for (unsigned int j = 0; j < sublayers[i]; ++j)
_planes.push_back(_planes.back() + layer_size);
}
}
if (_planes.size() < 2)
paramError("plane_coordinates", "Size of 'plane_coordinates' should be at least two");
}

std::unique_ptr<MeshBase>
PlaneIDMeshGenerator::generate()
{
std::unique_ptr<MeshBase> mesh = std::move(_input);
if (mesh->mesh_dimension() < _axis_index + 1)
paramError("plane_axis", "PlaneIDMeshGenerator must operate on a proper plane axis");
unsigned int extra_id_index = 0;
if (!mesh->has_elem_integer(_element_id_name))
extra_id_index = mesh->add_elem_integer(_element_id_name);
else
{
extra_id_index = mesh->get_elem_integer_index(_element_id_name);
paramWarning(
"id_name", "An element integer with the name '", _element_id_name, "' already exists");
}

const Real tol = getParam<Real>("tolerance");
for (auto & elem : mesh->active_element_ptr_range())
{
const int layer_id = getPlaneID(elem->centroid());
for (unsigned int i = 0; i < elem->n_nodes(); ++i)
{
const Point & p = elem->point(i) - (elem->point(i) - elem->centroid()) * tol;
if (getPlaneID(p) != layer_id)
mooseError("Element at ", elem->centroid(), " is cut by the planes");
}
elem->set_extra_integer(extra_id_index, layer_id);
}
return dynamic_pointer_cast<MeshBase>(mesh);
}

int
PlaneIDMeshGenerator::getPlaneID(const Point & p) const
{
// check the input point is located within the plane range
if (p(_axis_index) < _planes[0] || p(_axis_index) > _planes.back())
mooseError("The planes do not cover element at ", p);
int id = 0;
// looping over each plane to find a plane ID corresponding to the input point
for (unsigned int layer_id = 0; layer_id < _planes.size() - 1; ++layer_id)
{
// check the input point is located between _planes[layer_id] and _planes[layer_id + 1]
if (_planes[layer_id] < p(_axis_index) && _planes[layer_id + 1] >= p(_axis_index))
{
id = layer_id;
break;
}
}
return id;
}
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
@@ -0,0 +1,48 @@
[Mesh]
[grid]
type = CartesianMeshGenerator
dim = 2
dx = '5.0 10.0 '
ix = '1 4'
dy = '5.0 10.0 '
iy = '1 4'
[]

[plane_id_gen]
type = PlaneIDMeshGenerator
input = 'grid'
plane_coordinates = '0.0 5.0 15.0'
num_ids_per_plane = ' 1 2'
plane_axis = 'x'
id_name = 'plane_id'
[]
[]


[Executioner]
type = Steady
[]

[Problem]
solve = false
[]

[AuxVariables]
[plane_id]
family = MONOMIAL
order = CONSTANT
[]
[]

[AuxKernels]
[set_plane_id]
type = ExtraElementIDAux
variable = plane_id
extra_id_name = plane_id
[]
[]

[Outputs]
exodus = true
execute_on = timestep_end
[]
@@ -0,0 +1,56 @@
[Mesh]
[pin2d]
type = PolygonConcentricCircleMeshGenerator
preserve_volumes = true
ring_radii = 0.4
ring_intervals = 2
background_intervals = 1
num_sides = 6
num_sectors_per_side = '2 2 2 2 2 2'
polygon_size = 0.5
[]
[pin3d]
type = FancyExtruderGenerator
input = 'pin2d'
heights = '5.0 5.0 5.0'
direction = '0 0 1'
num_layers = '2 2 2'
[]
[pin3d_id]
type = PlaneIDMeshGenerator
input = 'pin3d'
plane_coordinates = '0.0 5.0 10.0 15.0'
num_ids_per_plane = ' 1 2 1'
plane_axis = 'z'
id_name = 'plane_id'
[]
[]


[Executioner]
type = Steady
[]

[Problem]
solve = false
[]

[AuxVariables]
[plane_id]
family = MONOMIAL
order = CONSTANT
[]
[]

[AuxKernels]
[set_plane_id]
type = ExtraElementIDAux
variable = plane_id
extra_id_name = plane_id
[]
[]

[Outputs]
exodus = true
execute_on = timestep_end
[]
@@ -0,0 +1,28 @@
[Tests]
[2d_cartesian_grid]
requirement = 'The system shall support the generation of plane IDs for 2D Cartesian grid'
[x_dir]
type = 'Exodiff'
input = 'plane_id_grid2d.i'
exodiff = 'plane_id_grid2d_x_dir_out.e'
cli_args = "Outputs/file_base='plane_id_grid2d_x_dir_out'"
detail = 'by taking x-direction as the plane axis'
recover = false
[]
[y_dir]
type = 'Exodiff'
input = 'plane_id_grid2d.i'
exodiff = 'plane_id_grid2d_y_dir_out.e'
cli_args = "Mesh/plane_id_gen/plane_axis='y' Outputs/file_base='plane_id_grid2d_y_dir_out'"
detail = 'by taking y-direction as the plane axis'
recover = false
[]
[]
[3d_extruded_mesh]
type = 'Exodiff'
input = 'plane_id_pin3d.i'
exodiff = 'plane_id_pin3d_out.e'
requirement = 'The system shall support the generation of plane IDs for 3D extruded mesh'
recover = false
[]
[]

0 comments on commit 0a58ec5

Please sign in to comment.