Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
3 : //*
4 : //* All rights reserved, see COPYRIGHT for full restrictions
5 : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 : //*
7 : //* Licensed under LGPL 2.1, please see LICENSE for details
8 : //* https://www.gnu.org/licenses/lgpl-2.1.html
9 :
10 : #include "HexagonalGridDivision.h"
11 : #include "MooseMesh.h"
12 : #include "HexagonalLatticeUtils.h"
13 : #include "Positions.h"
14 :
15 : #include "libmesh/elem.h"
16 :
17 : registerMooseObject("ReactorApp", HexagonalGridDivision);
18 :
19 : InputParameters
20 62 : HexagonalGridDivision::validParams()
21 : {
22 62 : InputParameters params = MeshDivision::validParams();
23 62 : params.addClassDescription(
24 : "Divide the mesh along a hexagonal grid. Numbering of pin divisions increases first "
25 : "counterclockwise, then expanding outwards from the inner ring, then axially. "
26 : "Inner-numbering is within a radial ring, outer-numbering is axial divisions");
27 :
28 124 : params.addParam<Point>("center", "Center of the hexagonal grid");
29 124 : params.addParam<PositionsName>("center_positions", "Centers of the hexagonal grids");
30 :
31 124 : params.addRequiredRangeCheckedParam<Real>(
32 : "lattice_flat_to_flat",
33 : "lattice_flat_to_flat>0",
34 : "Distance between two (inner) opposite sides of a lattice. Also known as bundle pitch or "
35 : "inner flat-to-flat distance");
36 124 : params.addRequiredRangeCheckedParam<Real>("pin_pitch", "pin_pitch>0", "Distance between pins");
37 :
38 124 : params.addRequiredParam<Real>("z_min", "Minimal axial extent of the lattice");
39 124 : params.addRequiredParam<Real>("z_max", "Maximum axial extent of the lattice");
40 124 : params.addRequiredRangeCheckedParam<unsigned int>("nr", "nr>0", "Number of hexagonal rings");
41 124 : params.addRequiredRangeCheckedParam<unsigned int>("nz", "nz>0", "Number of divisions in Z");
42 124 : params.addParam<bool>(
43 : "assign_domain_outside_grid_to_border",
44 124 : false,
45 : "Whether to map the domain outside the grid back to the border of the grid");
46 124 : params.addParam<Real>("rotation_around_axis",
47 124 : 0.,
48 : "Rotation angle to apply to the underlying hexagonal lattice (in degrees)");
49 :
50 62 : return params;
51 0 : }
52 :
53 34 : HexagonalGridDivision::HexagonalGridDivision(const InputParameters & parameters)
54 : : MeshDivision(parameters),
55 90 : _center(isParamValid("center") ? getParam<Point>("center") : Point(0, 0, 0)),
56 34 : _center_positions(
57 34 : isParamValid("center_positions")
58 58 : ? &_fe_problem->getPositionsObject(getParam<PositionsName>("center_positions"))
59 : : nullptr),
60 68 : _lattice_flat_to_flat(getParam<Real>("lattice_flat_to_flat")),
61 68 : _pin_pitch(getParam<Real>("pin_pitch")),
62 68 : _z_axis_index(MooseEnum("X Y Z", "Z")),
63 68 : _min_z(getParam<Real>("z_min")),
64 68 : _max_z(getParam<Real>("z_max")),
65 68 : _nr(getParam<unsigned int>("nr")),
66 68 : _nz(getParam<unsigned int>("nz")),
67 102 : _outside_grid_counts_as_border(getParam<bool>("assign_domain_outside_grid_to_border"))
68 : {
69 34 : HexagonalGridDivision::initialize();
70 :
71 102 : if (!isParamValid("center") && !_center_positions)
72 0 : paramError("center", "A center must be provided, or a Positions object for the centers");
73 34 : if (_pin_pitch > _lattice_flat_to_flat)
74 0 : mooseError("lattice_flat_to_flat", "Pin pitch should be smaller than bundle pitch");
75 34 : if ((_nz > 1) && MooseUtils::absoluteFuzzyEqual(_max_z, _min_z))
76 0 : paramError("nz", "Subdivision number must be 1 if width is 0 in Z direction");
77 34 : }
78 :
79 : void
80 34 : HexagonalGridDivision::initialize()
81 : {
82 : // We make very large pins so they cover the entire position
83 68 : _hex_latt = std::make_unique<HexagonalLatticeUtils>(_lattice_flat_to_flat,
84 : _pin_pitch,
85 34 : _pin_pitch,
86 68 : 0.,
87 34 : 1.,
88 34 : _nr,
89 34 : _z_axis_index,
90 : getParam<Real>("rotation_around_axis"));
91 :
92 34 : if (!_center_positions)
93 22 : setNumDivisions(_hex_latt->totalPins(_nr) * _nz);
94 : else
95 12 : setNumDivisions(_center_positions->getNumPositions() * _hex_latt->totalPins(_nr) * _nz);
96 :
97 : // Check that the grid is well-defined
98 34 : if (_center_positions)
99 : {
100 12 : const Real min_center_dist = _center_positions->getMinDistanceBetweenPositions();
101 : // Note that if the positions are not aligned on a hexagonal lattice themselves,
102 : // this bound is not sufficiently strict. The simplest example would be non-coplanar
103 : // points, which can be a great distance away axially but be on the same axis
104 12 : if (MooseUtils::absoluteFuzzyGreaterThan(_lattice_flat_to_flat, min_center_dist))
105 0 : mooseWarning(
106 : "Hexagonal grids centered on the positions are too close to each other (min distance: ",
107 : min_center_dist,
108 : "), closer than the extent of each grid (",
109 : _lattice_flat_to_flat,
110 : "). Mesh division is ill-defined ");
111 : }
112 34 : }
113 :
114 : unsigned int
115 28008 : HexagonalGridDivision::divisionIndex(const Elem & elem) const
116 : {
117 28008 : return divisionIndex(elem.vertex_average());
118 : }
119 :
120 : unsigned int
121 28008 : HexagonalGridDivision::divisionIndex(const Point & pt) const
122 : {
123 : unsigned int offset = 0;
124 :
125 : // Get point in the coordinates of the lattice. This can involve a projection due to
126 : // the axis of the lattice, or simply a translation if there are lattices distributed
127 : // using positions
128 : Point pc;
129 28008 : if (_center_positions)
130 : {
131 : // If dividing using positions, find the closest position and
132 : // look at the relative position of the point compared to that position
133 8568 : const bool initial = _fe_problem->getCurrentExecuteOnFlag() == EXEC_INITIAL;
134 8568 : const auto nearest_grid_center_index = _center_positions->getNearestPositionIndex(pt, initial);
135 8568 : offset = nearest_grid_center_index * _hex_latt->totalPins(_nr) * _nz;
136 : const auto nearest_grid_center =
137 8568 : _center_positions->getPosition(nearest_grid_center_index, initial);
138 :
139 : // Project in local hexagonal grid
140 8568 : pc = pt - nearest_grid_center;
141 : }
142 : else
143 19440 : pc = pt - _center;
144 :
145 : // Get radial division index, using the channel as the pins are 0-radius
146 : // The logic in get pin index requires getting the point in the plane of the pin centers
147 28008 : auto ir = _hex_latt->pinIndex(pc);
148 : const auto n_pins = _hex_latt->nPins();
149 :
150 28008 : if (!_outside_grid_counts_as_border)
151 : {
152 28008 : if (ir == n_pins)
153 8532 : return MooseMeshDivision::INVALID_DIVISION_INDEX;
154 19476 : if (MooseUtils::absoluteFuzzyLessThan(pc(_z_axis_index), _min_z) ||
155 : MooseUtils::absoluteFuzzyGreaterThan(pc(_z_axis_index), _max_z))
156 0 : return MooseMeshDivision::INVALID_DIVISION_INDEX;
157 : }
158 :
159 : // If too far from the grid to have a valid radial index, use the closest pin
160 19476 : if (ir == n_pins)
161 0 : ir = _hex_latt->closestPinIndex(pc);
162 :
163 : // Find axial index
164 19476 : const auto not_found = MooseMeshDivision::INVALID_DIVISION_INDEX;
165 : auto iz = not_found;
166 35001 : for (const auto jz : make_range(_nz + 1))
167 : {
168 35001 : const auto border_z = _min_z + (_max_z - _min_z) * jz / _nz;
169 35001 : if (jz > 0 && jz < _nz && MooseUtils::absoluteFuzzyEqual(border_z, pc(_z_axis_index)))
170 0 : mooseWarning(
171 0 : "Querying the division index for a point of a boundary between two regions in Z: " +
172 0 : Moose::stringify(pt),
173 : ", in local hex grid frame: ",
174 0 : Moose::stringify(pc));
175 35001 : if (border_z >= pc(_z_axis_index))
176 : {
177 19476 : iz = (jz > 0) ? jz - 1 : 0;
178 : break;
179 : }
180 : }
181 :
182 : // Look on the top of the grid
183 19476 : if (MooseUtils::absoluteFuzzyGreaterEqual(pc(_z_axis_index), _max_z))
184 7056 : iz = _nz - 1;
185 :
186 : // Handle edge case on widths
187 19476 : if (iz == not_found && MooseUtils::absoluteFuzzyEqual(_max_z - _min_z, 0))
188 : iz = 0;
189 : mooseAssert(ir != not_found, "We should have found a mesh division bin radially");
190 : mooseAssert(iz != not_found, "We should have found a mesh division bin in Z");
191 :
192 19476 : const auto n_radial = _hex_latt->totalPins(_nr);
193 19476 : return offset + ir + iz * n_radial;
194 : }
|