https://mooseframework.inl.gov
Public Types | Public Member Functions | Static Public Member Functions | Public Attributes | Static Public Attributes | Protected Member Functions | Static Protected Member Functions | Protected Attributes | List of all members
SCMDetailedQuadSubChannelMeshGenerator Class Reference

Mesh generator that builds a 3D mesh representing quadrilateral subchannels. More...

#include <SCMDetailedQuadSubChannelMeshGenerator.h>

Inheritance diagram for SCMDetailedQuadSubChannelMeshGenerator:
[legend]

Public Types

typedef DataFileName DataFileParameterType
 

Public Member Functions

 SCMDetailedQuadSubChannelMeshGenerator (const InputParameters &parameters)
 
virtual std::unique_ptr< MeshBase > generate () override
 
std::unique_ptr< MeshBase > generateInternal ()
 
const std::set< MeshGeneratorName > & getRequestedMeshGenerators () const
 
const std::set< MeshGeneratorName > & getRequestedMeshGeneratorsForSub () const
 
void addParentMeshGenerator (const MeshGenerator &mg, const AddParentChildKey)
 
void addChildMeshGenerator (const MeshGenerator &mg, const AddParentChildKey)
 
const std::set< const MeshGenerator *, Comparator > & getParentMeshGenerators () const
 
const std::set< const MeshGenerator *, Comparator > & getChildMeshGenerators () const
 
const std::set< const MeshGenerator *, Comparator > & getSubMeshGenerators () const
 
bool isParentMeshGenerator (const MeshGeneratorName &name, const bool direct=true) const
 
bool isChildMeshGenerator (const MeshGeneratorName &name, const bool direct=true) const
 
bool isNullMeshName (const MeshGeneratorName &name) const
 
bool hasSaveMesh () const
 
bool hasOutput () const
 
const std::string & getSavedMeshName () const
 
bool hasGenerateData () const
 
bool isDataOnly () const
 
virtual bool enabled () const
 
std::shared_ptr< MooseObjectgetSharedPtr ()
 
std::shared_ptr< const MooseObjectgetSharedPtr () const
 
MooseAppgetMooseApp () const
 
const std::string & type () const
 
const std::string & name () const
 
std::string typeAndName () const
 
MooseObjectParameterName uniqueParameterName (const std::string &parameter_name) const
 
MooseObjectName uniqueName () const
 
const InputParametersparameters () const
 
const hit::Node * getHitNode () const
 
bool hasBase () const
 
const std::string & getBase () const
 
const T & getParam (const std::string &name) const
 
std::vector< std::pair< T1, T2 > > getParam (const std::string &param1, const std::string &param2) const
 
const T * queryParam (const std::string &name) const
 
const T & getRenamedParam (const std::string &old_name, const std::string &new_name) const
 
getCheckedPointerParam (const std::string &name, const std::string &error_string="") const
 
bool isParamValid (const std::string &name) const
 
bool isParamSetByUser (const std::string &name) const
 
void connectControllableParams (const std::string &parameter, const std::string &object_type, const std::string &object_name, const std::string &object_parameter) const
 
void paramError (const std::string &param, Args... args) const
 
void paramWarning (const std::string &param, Args... args) const
 
void paramInfo (const std::string &param, Args... args) const
 
std::string messagePrefix (const bool hit_prefix=true) const
 
std::string errorPrefix (const std::string &) const
 
void mooseError (Args &&... args) const
 
void mooseDocumentedError (const std::string &repo_name, const unsigned int issue_num, Args &&... args) const
 
void mooseErrorNonPrefixed (Args &&... args) const
 
void mooseWarning (Args &&... args) const
 
void mooseWarningNonPrefixed (Args &&... args) const
 
void mooseDeprecated (Args &&... args) const
 
void mooseInfo (Args &&... args) const
 
void callMooseError (std::string msg, const bool with_prefix, const hit::Node *node=nullptr) const
 
std::string getDataFileName (const std::string &param) const
 
std::string getDataFileNameByName (const std::string &relative_path) const
 
std::string getDataFilePath (const std::string &relative_path) const
 
const Parallel::Communicator & comm () const
 
processor_id_type n_processors () const
 
processor_id_type processor_id () const
 

Static Public Member Functions

static InputParameters validParams ()
 
static bool hasGenerateData (const InputParameters &params)
 
static void callMooseError (MooseApp *const app, const InputParameters &params, std::string msg, const bool with_prefix, const hit::Node *node)
 
static void setHasGenerateData (InputParameters &params)
 

Public Attributes

const ConsoleStream _console
 

Static Public Attributes

static const std::string data_only_param
 
static const std::string type_param
 
static const std::string name_param
 
static const std::string unique_name_param
 
static const std::string app_param
 
static const std::string moose_base_param
 
static constexpr auto SYSTEM
 
static constexpr auto NAME
 

Protected Member Functions

EChannelType getSubchannelType (unsigned int index) const
 
std::vector< RealgetSubchannelPosition (unsigned int i)
 
virtual void generateData ()
 
T & copyMeshProperty (const std::string &target_data_name, const std::string &source_data_name, const std::string &source_mesh)
 
T & copyMeshProperty (const std::string &source_data_name, const std::string &source_mesh)
 
std::unique_ptr< MeshBase > & getMesh (const std::string &param_name, const bool allow_invalid=false)
 
std::vector< std::unique_ptr< MeshBase > *> getMeshes (const std::string &param_name)
 
std::unique_ptr< MeshBase > & getMeshByName (const MeshGeneratorName &mesh_generator_name)
 
std::vector< std::unique_ptr< MeshBase > *> getMeshesByName (const std::vector< MeshGeneratorName > &mesh_generator_names)
 
void declareMeshForSub (const std::string &param_name)
 
void declareMeshesForSub (const std::string &param_name)
 
void declareMeshForSubByName (const MeshGeneratorName &mesh_generator_name)
 
void declareMeshesForSubByName (const std::vector< MeshGeneratorName > &mesh_generator_names)
 
std::unique_ptr< MeshBase > buildMeshBaseObject (unsigned int dim=libMesh::invalid_uint)
 
std::unique_ptr< ReplicatedMesh > buildReplicatedMesh (unsigned int dim=libMesh::invalid_uint)
 
std::unique_ptr< DistributedMesh > buildDistributedMesh (unsigned int dim=libMesh::invalid_uint)
 
void addMeshSubgenerator (const std::string &type, const std::string &name, Ts... extra_input_parameters)
 
void addMeshSubgenerator (const std::string &type, const std::string &name, InputParameters params)
 
void declareNullMeshName (const MeshGeneratorName &name)
 
const T & getMeshProperty (const std::string &data_name, const std::string &prefix)
 
const T & getMeshProperty (const std::string &data_name)
 
bool hasMeshProperty (const std::string &data_name, const std::string &prefix) const
 
bool hasMeshProperty (const std::string &data_name, const std::string &prefix) const
 
bool hasMeshProperty (const std::string &data_name) const
 
bool hasMeshProperty (const std::string &data_name) const
 
std::string meshPropertyName (const std::string &data_name) const
 
T & declareMeshProperty (const std::string &data_name, Args &&... args)
 
T & declareMeshProperty (const std::string &data_name, const T &data_value)
 
T & declareMeshProperty (const std::string &data_name, Args &&... args)
 
T & declareMeshProperty (const std::string &data_name, const T &data_value)
 
T & setMeshProperty (const std::string &data_name, Args &&... args)
 
T & setMeshProperty (const std::string &data_name, const T &data_value)
 
T & setMeshProperty (const std::string &data_name, Args &&... args)
 
T & setMeshProperty (const std::string &data_name, const T &data_value)
 

Static Protected Member Functions

static std::string meshPropertyName (const std::string &data_name, const std::string &prefix)
 

Protected Attributes

const Real _unheated_length_entry
 unheated length of the fuel Pin at the entry of the assembly More...
 
const Real _heated_length
 heated length of the fuel Pin More...
 
const Real _unheated_length_exit
 unheated length of the fuel Pin at the exit of the assembly More...
 
std::vector< Real_z_grid
 axial location of nodes More...
 
const Real _pitch
 Distance between the neighbor fuel pins, pitch. More...
 
const Real _pin_diameter
 fuel Pin diameter More...
 
const unsigned int _n_cells
 Number of cells in the axial direction. More...
 
const unsigned int _nx
 Number of subchannels in the x direction. More...
 
const unsigned int _ny
 Number of subchannels in the y direction. More...
 
unsigned int _n_channels
 Total number of subchannels. More...
 
const Real _side_gap
 The side gap, not to be confused with the gap between pins, this refers to the gap next to the duct or else the distance between the subchannel centroid to the duct wall. More...
 
std::vector< EChannelType_subch_type
 Subchannel type. More...
 
std::vector< std::vector< Real > > _subchannel_position
 x,y coordinates of the subchannel centroids More...
 
const unsigned int_block_id
 Subdomain ID used for the mesh block. More...
 
MooseMesh *const _mesh
 
const bool & _enabled
 
MooseApp_app
 
Factory_factory
 
ActionFactory_action_factory
 
const std::string & _type
 
const std::string & _name
 
const InputParameters_pars
 
const Parallel::Communicator & _communicator
 

Detailed Description

Mesh generator that builds a 3D mesh representing quadrilateral subchannels.

Definition at line 18 of file SCMDetailedQuadSubChannelMeshGenerator.h.

Constructor & Destructor Documentation

◆ SCMDetailedQuadSubChannelMeshGenerator()

SCMDetailedQuadSubChannelMeshGenerator::SCMDetailedQuadSubChannelMeshGenerator ( const InputParameters parameters)

Definition at line 46 of file SCMDetailedQuadSubChannelMeshGenerator.C.

49  _unheated_length_entry(getParam<Real>("unheated_length_entry")),
50  _heated_length(getParam<Real>("heated_length")),
51  _unheated_length_exit(getParam<Real>("unheated_length_exit")),
52  _pitch(getParam<Real>("pitch")),
53  _pin_diameter(getParam<Real>("pin_diameter")),
54  _n_cells(getParam<unsigned int>("n_cells")),
55  _nx(getParam<unsigned int>("nx")),
56  _ny(getParam<unsigned int>("ny")),
57  _n_channels(_nx * _ny),
58  _side_gap(getParam<Real>("side_gap")),
59  _block_id(getParam<unsigned int>("block_id"))
60 {
62  Real dz = L / _n_cells;
63  for (unsigned int i = 0; i < _n_cells + 1; i++)
64  _z_grid.push_back(dz * i);
65 
67  for (unsigned int i = 0; i < _n_channels; i++)
68  {
69  _subchannel_position[i].reserve(3);
70  for (unsigned int j = 0; j < 3; j++)
71  {
72  _subchannel_position.at(i).push_back(0.0);
73  }
74  }
75 
76  _subch_type.resize(_n_channels);
77  for (unsigned int iy = 0; iy < _ny; iy++)
78  {
79  for (unsigned int ix = 0; ix < _nx; ix++)
80  {
81  unsigned int i_ch = _nx * iy + ix;
82  bool is_corner = (ix == 0 && iy == 0) || (ix == _nx - 1 && iy == 0) ||
83  (ix == 0 && iy == _ny - 1) || (ix == _nx - 1 && iy == _ny - 1);
84  bool is_edge = (ix == 0 || iy == 0 || ix == _nx - 1 || iy == _ny - 1);
85 
86  if (is_corner)
88  else if (is_edge)
90  else
92 
93  // set the subchannel positions so that the center of the assembly is the zero point
94  Real offset_x = (_nx - 1) * _pitch / 2.0;
95  Real offset_y = (_ny - 1) * _pitch / 2.0;
96  _subchannel_position[i_ch][0] = _pitch * ix - offset_x;
97  _subchannel_position[i_ch][1] = _pitch * iy - offset_y;
98  }
99  }
100 }
const Real _pitch
Distance between the neighbor fuel pins, pitch.
const unsigned int _n_cells
Number of cells in the axial direction.
const Real _heated_length
heated length of the fuel Pin
std::vector< Real > _z_grid
axial location of nodes
unsigned int _n_channels
Total number of subchannels.
const InputParameters & parameters() const
MeshGenerator(const InputParameters &parameters)
const unsigned int _nx
Number of subchannels in the x direction.
std::vector< std::vector< Real > > _subchannel_position
x,y coordinates of the subchannel centroids
const Real _side_gap
The side gap, not to be confused with the gap between pins, this refers to the gap next to the duct o...
const unsigned int _ny
Number of subchannels in the y direction.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const Real _unheated_length_exit
unheated length of the fuel Pin at the exit of the assembly
std::vector< EChannelType > _subch_type
Subchannel type.
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
const Real _unheated_length_entry
unheated length of the fuel Pin at the entry of the assembly
const unsigned int & _block_id
Subdomain ID used for the mesh block.

Member Function Documentation

◆ generate()

std::unique_ptr< MeshBase > SCMDetailedQuadSubChannelMeshGenerator::generate ( )
overridevirtual

Implements MeshGenerator.

Definition at line 103 of file SCMDetailedQuadSubChannelMeshGenerator.C.

104 {
105  auto mesh_base = buildMeshBaseObject();
106  BoundaryInfo & boundary_info = mesh_base->get_boundary_info();
107  mesh_base->set_spatial_dimension(3);
108  // Define the resolution (the number of points used to represent a circle).
109  // This must be divisible by 4.
110  const unsigned int theta_res = 16; // TODO: parameterize
111  // Compute the number of points needed to represent one quarter of a circle.
112  const unsigned int points_per_quad = theta_res / 4 + 1;
113 
114  // Compute the points needed to represent one axial cross-flow of a subchannel.
115  // For the center subchannel (sc) there is one center point plus the points from 4 intersecting
116  // circles.
117  const unsigned int points_per_center = points_per_quad * 4 + 1;
118  // For the corner sc there is one center point plus the points from 1 intersecting circle plus 3
119  // corners
120  const unsigned int points_per_corner = points_per_quad * 1 + 1 + 3;
121  // For the side sc there is one center point plus the points from 2 intersecting circles plus 2
122  // corners
123  const unsigned int points_per_side = points_per_quad * 2 + 1 + 2;
124 
125  // Compute the number of elements (Prism6) which combined base creates the sub-channel
126  // cross-section
127  const unsigned int elems_per_center = theta_res + 4;
128  const unsigned int elems_per_corner = theta_res / 4 + 4;
129  const unsigned int elems_per_side = theta_res / 2 + 4;
130 
131  // specify number and type of sub-channel
132  unsigned int n_center, n_side, n_corner;
133  if (_n_channels == 2)
134  {
135  n_corner = 0;
136  n_side = _n_channels;
137  n_center = _n_channels - n_side - n_corner;
138  }
139  else if (_n_channels > 2 && (_ny == 1 || _nx == 1))
140  {
141  n_corner = 0;
142  n_side = 2;
143  n_center = _n_channels - n_side - n_corner;
144  }
145  else
146  {
147  n_corner = 4;
148  n_side = 2 * (_ny - 2) + 2 * (_nx - 2);
149  n_center = _n_channels - n_side - n_corner;
150  }
151 
152  // Compute the total number of points and elements.
153  const unsigned int points_per_level =
154  n_corner * points_per_corner + n_side * points_per_side + n_center * points_per_center;
155  const unsigned int elems_per_level =
156  n_corner * elems_per_corner + n_side * elems_per_side + n_center * elems_per_center;
157  const unsigned int n_points = points_per_level * (_n_cells + 1);
158  const unsigned int n_elems = elems_per_level * _n_cells;
159  mesh_base->reserve_nodes(n_points);
160  mesh_base->reserve_elem(n_elems);
161  // Build an array of points arranged in a circle on the xy-plane. (last and first node overlap)
162  const Real radius = _pin_diameter / 2.0;
163  std::array<Point, theta_res + 1> circle_points;
164  {
165  Real theta = 0;
166  for (unsigned int i = 0; i < theta_res + 1; i++)
167  {
168  circle_points[i](0) = radius * std::cos(theta);
169  circle_points[i](1) = radius * std::sin(theta);
170  theta += 2 * M_PI / theta_res;
171  }
172  }
173  // Define "quadrant center" reference points. These will be the centers of
174  // the 4 circles that represent the fuel pins. These centers are
175  // offset a little bit so that in the final mesh, there is a tiny gap between
176  // neighboring subchannel cells. That allows us to easily map a solution to
177  // this detailed mesh with a nearest-neighbor search.
178  const Real shrink_factor = 0.99999;
179  std::array<Point, 4> quadrant_centers;
180  quadrant_centers[0] = Point(_pitch * 0.5 * shrink_factor, _pitch * 0.5 * shrink_factor, 0);
181  quadrant_centers[1] = Point(-_pitch * 0.5 * shrink_factor, _pitch * 0.5 * shrink_factor, 0);
182  quadrant_centers[2] = Point(-_pitch * 0.5 * shrink_factor, -_pitch * 0.5 * shrink_factor, 0);
183  quadrant_centers[3] = Point(_pitch * 0.5 * shrink_factor, -_pitch * 0.5 * shrink_factor, 0);
184 
185  const unsigned int m = theta_res / 4;
186  // Build an array of points that represent a cross section of a center subchannel
187  // cell. The points are ordered in this fashion:
188  // 4 3
189  // 6 5 2 1
190  // 0
191  // 7 8 * *
192  // 9 *
193  std::array<Point, points_per_center> center_points;
194  {
195  unsigned int start;
196  for (unsigned int i = 0; i < 4; i++)
197  {
198  if (i == 0)
199  start = 3 * m;
200  if (i == 1)
201  start = 4 * m;
202  if (i == 2)
203  start = 1 * m;
204  if (i == 3)
205  start = 2 * m;
206  for (unsigned int ii = 0; ii < points_per_quad; ii++)
207  {
208  auto c_pt = circle_points[start - ii];
209  center_points[i * points_per_quad + ii + 1] = quadrant_centers[i] + c_pt;
210  }
211  }
212  }
213 
214  // Build an array of points that represent a cross section of a top left corner subchannel
215  // cell. The points are ordered in this fashion:
216  // 5 4
217  //
218  // 0
219  // 2 3
220  // 6 1
221  std::array<Point, points_per_corner> tl_corner_points;
222  {
223  for (unsigned int ii = 0; ii < points_per_quad; ii++)
224  {
225  auto c_pt = circle_points[2 * m - ii];
226  tl_corner_points[ii + 1] = quadrant_centers[3] + c_pt;
227  }
228  tl_corner_points[points_per_quad + 1] = Point(_pitch * 0.5 * shrink_factor, _side_gap, 0);
229  tl_corner_points[points_per_quad + 2] = Point(-_side_gap, _side_gap, 0);
230  tl_corner_points[points_per_quad + 3] = Point(-_side_gap, -_pitch * 0.5 * shrink_factor, 0);
231  }
232 
233  // Build an array of points that represent a cross section of a top right corner subchannel
234  // cell. The points are ordered in this fashion:
235  // 6 5
236  //
237  // 0
238  // 1 2
239  // 3 4
240  std::array<Point, points_per_corner> tr_corner_points;
241  {
242  for (unsigned int ii = 0; ii < points_per_quad; ii++)
243  {
244  auto c_pt = circle_points[m - ii];
245  tr_corner_points[ii + 1] = quadrant_centers[2] + c_pt;
246  }
247  tr_corner_points[points_per_quad + 1] = Point(_side_gap, -_pitch * 0.5 * shrink_factor, 0);
248  tr_corner_points[points_per_quad + 2] = Point(_side_gap, _side_gap, 0);
249  tr_corner_points[points_per_quad + 3] = Point(-_pitch * 0.5 * shrink_factor, _side_gap, 0);
250  }
251 
252  // Build an array of points that represent a cross section of a bottom left corner subchannel
253  // cell. The points are ordered in this fashion:
254  // 4 3
255  // 2 1
256  // 0
257  //
258  // 5 6
259  std::array<Point, points_per_corner> bl_corner_points;
260  {
261  for (unsigned int ii = 0; ii < points_per_quad; ii++)
262  {
263  auto c_pt = circle_points[3 * m - ii];
264  bl_corner_points[ii + 1] = quadrant_centers[0] + c_pt;
265  }
266  bl_corner_points[points_per_quad + 1] = Point(-_side_gap, _pitch * 0.5 * shrink_factor, 0);
267  bl_corner_points[points_per_quad + 2] = Point(-_side_gap, -_side_gap, 0);
268  bl_corner_points[points_per_quad + 3] = Point(_pitch * 0.5 * shrink_factor, -_side_gap, 0);
269  }
270 
271  // Build an array of points that represent a cross section of a bottom right corner subchannel
272  // cell. The points are ordered in this fashion:
273  // 1 6
274  // 3 2
275  // 0
276  //
277  // 4 5
278  std::array<Point, points_per_corner> br_corner_points;
279  {
280  for (unsigned int ii = 0; ii < points_per_quad; ii++)
281  {
282  auto c_pt = circle_points[4 * m - ii];
283  br_corner_points[ii + 1] = quadrant_centers[1] + c_pt;
284  }
285  br_corner_points[points_per_quad + 1] = Point(-_pitch * 0.5 * shrink_factor, -_side_gap, 0);
286  br_corner_points[points_per_quad + 2] = Point(_side_gap, -_side_gap, 0);
287  br_corner_points[points_per_quad + 3] = Point(_side_gap, _pitch * 0.5 * shrink_factor, 0);
288  }
289 
290  // Build an array of points that represent a cross section of a top side subchannel
291  // cell. The points are ordered in this fashion:
292  // 8 7
293  //
294  // 0
295  // 1 2 5 6
296  // 3 4
297  std::array<Point, points_per_side> top_points;
298  {
299  for (unsigned int ii = 0; ii < points_per_quad; ii++)
300  {
301  auto c_pt = circle_points[m - ii];
302  top_points[ii + 1] = quadrant_centers[2] + c_pt;
303  }
304  for (unsigned int ii = 0; ii < points_per_quad; ii++)
305  {
306  auto c_pt = circle_points[2 * m - ii];
307  top_points[points_per_quad + ii + 1] = quadrant_centers[3] + c_pt;
308  }
309  top_points[2 * points_per_quad + 1] = Point(_pitch * 0.5 * shrink_factor, _side_gap, 0);
310  top_points[2 * points_per_quad + 2] = Point(-_pitch * 0.5 * shrink_factor, _side_gap, 0);
311  }
312 
313  // Build an array of points that represent a cross section of a left side subchannel
314  // cell. The points are ordered in this fashion:
315  // 7 6
316  // 5 4
317  // 0
318  // 2 3
319  // 8 1
320  std::array<Point, points_per_side> left_points;
321  {
322  for (unsigned int ii = 0; ii < points_per_quad; ii++)
323  {
324  auto c_pt = circle_points[2 * m - ii];
325  left_points[ii + 1] = quadrant_centers[3] + c_pt;
326  }
327  for (unsigned int ii = 0; ii < points_per_quad; ii++)
328  {
329  auto c_pt = circle_points[3 * m - ii];
330  left_points[points_per_quad + ii + 1] = quadrant_centers[0] + c_pt;
331  }
332  left_points[2 * points_per_quad + 1] = Point(-_side_gap, _pitch * 0.5 * shrink_factor, 0);
333  left_points[2 * points_per_quad + 2] = Point(-_side_gap, -_pitch * 0.5 * shrink_factor, 0);
334  }
335 
336  // Build an array of points that represent a cross section of a bottom side subchannel
337  // cell. The points are ordered in this fashion:
338  // 4 3
339  // 6 5 2 1
340  // 0
341  //
342  // 7 8
343  std::array<Point, points_per_side> bottom_points;
344  {
345  for (unsigned int ii = 0; ii < points_per_quad; ii++)
346  {
347  auto c_pt = circle_points[3 * m - ii];
348  bottom_points[ii + 1] = quadrant_centers[0] + c_pt;
349  }
350  for (unsigned int ii = 0; ii < points_per_quad; ii++)
351  {
352  auto c_pt = circle_points[4 * m - ii];
353  bottom_points[points_per_quad + ii + 1] = quadrant_centers[1] + c_pt;
354  }
355  bottom_points[2 * points_per_quad + 1] = Point(-_pitch * 0.5 * shrink_factor, -_side_gap, 0);
356  bottom_points[2 * points_per_quad + 2] = Point(_pitch * 0.5 * shrink_factor, -_side_gap, 0);
357  }
358 
359  // Build an array of points that represent a cross section of a right side subchannel
360  // cell. The points are ordered in this fashion:
361  // 1 8
362  // 3 2
363  // 0
364  // 4 5
365  // 6 7
366  std::array<Point, points_per_side> right_points;
367  {
368  for (unsigned int ii = 0; ii < points_per_quad; ii++)
369  {
370  auto c_pt = circle_points[4 * m - ii];
371  right_points[ii + 1] = quadrant_centers[1] + c_pt;
372  }
373  for (unsigned int ii = 0; ii < points_per_quad; ii++)
374  {
375  auto c_pt = circle_points[m - ii];
376  right_points[points_per_quad + ii + 1] = quadrant_centers[2] + c_pt;
377  }
378  right_points[2 * points_per_quad + 1] = Point(_side_gap, -_pitch * 0.5 * shrink_factor, 0);
379  right_points[2 * points_per_quad + 2] = Point(_side_gap, _pitch * 0.5 * shrink_factor, 0);
380  }
381 
382  // Add the points to the mesh.
383  if (_n_channels == 2)
384  {
385  unsigned int node_id = 0;
386  Real offset_x = (_nx - 1) * _pitch / 2.0;
387  Real offset_y = (_ny - 1) * _pitch / 2.0;
388  for (unsigned int iy = 0; iy < _ny; iy++)
389  {
390  Point y0 = {0, _pitch * iy - offset_y, 0};
391  for (unsigned int ix = 0; ix < _nx; ix++)
392  {
393  Point x0 = {_pitch * ix - offset_x, 0, 0};
394  for (auto z : _z_grid)
395  {
396  Point z0{0, 0, z};
397  if (_nx == 1 && iy == 0) // vertical orientation
398  {
399  for (unsigned int i = 0; i < points_per_side; i++)
400  mesh_base->add_point(bottom_points[i] + x0 + y0 + z0, node_id++);
401  }
402  if (_nx == 1 && iy == 1) // vertical orientation
403  {
404  for (unsigned int i = 0; i < points_per_side; i++)
405  mesh_base->add_point(top_points[i] + x0 + y0 + z0, node_id++);
406  }
407  if (_ny == 1 && ix == 0) // horizontal orientation
408  {
409  for (unsigned int i = 0; i < points_per_side; i++)
410  mesh_base->add_point(left_points[i] + x0 + y0 + z0, node_id++);
411  }
412  if (_ny == 1 && ix == 1) // horizontal orientation
413  {
414  for (unsigned int i = 0; i < points_per_side; i++)
415  mesh_base->add_point(right_points[i] + x0 + y0 + z0, node_id++);
416  }
417  }
418  }
419  }
420  }
421  else if (_n_channels > 2 && (_ny == 1 || _nx == 1))
422  {
423  unsigned int node_id = 0;
424  Real offset_x = (_nx - 1) * _pitch / 2.0;
425  Real offset_y = (_ny - 1) * _pitch / 2.0;
426  for (unsigned int iy = 0; iy < _ny; iy++)
427  {
428  Point y0 = {0, _pitch * iy - offset_y, 0};
429  for (unsigned int ix = 0; ix < _nx; ix++)
430  {
431  Point x0 = {_pitch * ix - offset_x, 0, 0};
432  for (auto z : _z_grid)
433  {
434  Point z0{0, 0, z};
435  if (_nx == 1)
436  {
437  if (iy == 0)
438  {
439  for (unsigned int i = 0; i < points_per_side; i++)
440  mesh_base->add_point(bottom_points[i] + x0 + y0 + z0, node_id++);
441  }
442  else if (iy == _ny - 1)
443  {
444  for (unsigned int i = 0; i < points_per_side; i++)
445  mesh_base->add_point(top_points[i] + x0 + y0 + z0, node_id++);
446  }
447  else
448  {
449  for (unsigned int i = 0; i < points_per_center; i++)
450  mesh_base->add_point(center_points[i] + x0 + y0 + z0, node_id++);
451  }
452  }
453  else if (_ny == 1)
454  {
455  if (ix == 0)
456  {
457  for (unsigned int i = 0; i < points_per_side; i++)
458  mesh_base->add_point(left_points[i] + x0 + y0 + z0, node_id++);
459  }
460  else if (ix == _nx - 1)
461  {
462  for (unsigned int i = 0; i < points_per_side; i++)
463  mesh_base->add_point(right_points[i] + x0 + y0 + z0, node_id++);
464  }
465  else
466  {
467  for (unsigned int i = 0; i < points_per_center; i++)
468  mesh_base->add_point(center_points[i] + x0 + y0 + z0, node_id++);
469  }
470  }
471  }
472  }
473  }
474  }
475  else
476  {
477  unsigned int node_id = 0;
478  Real offset_x = (_nx - 1) * _pitch / 2.0;
479  Real offset_y = (_ny - 1) * _pitch / 2.0;
480  for (unsigned int iy = 0; iy < _ny; iy++)
481  {
482  Point y0 = {0, _pitch * iy - offset_y, 0};
483  for (unsigned int ix = 0; ix < _nx; ix++)
484  {
485  Point x0 = {_pitch * ix - offset_x, 0, 0};
486  if (ix == 0 && iy == 0) // Bottom Left corner
487  {
488  for (auto z : _z_grid)
489  {
490  Point z0{0, 0, z};
491  for (unsigned int i = 0; i < points_per_corner; i++)
492  mesh_base->add_point(bl_corner_points[i] + x0 + y0 + z0, node_id++);
493  }
494  }
495  else if (ix == _nx - 1 && iy == 0) // Bottom right corner
496  {
497  for (auto z : _z_grid)
498  {
499  Point z0{0, 0, z};
500  for (unsigned int i = 0; i < points_per_corner; i++)
501  mesh_base->add_point(br_corner_points[i] + x0 + y0 + z0, node_id++);
502  }
503  }
504  else if (ix == 0 && iy == _ny - 1) // top Left corner
505  {
506  for (auto z : _z_grid)
507  {
508  Point z0{0, 0, z};
509  for (unsigned int i = 0; i < points_per_corner; i++)
510  mesh_base->add_point(tl_corner_points[i] + x0 + y0 + z0, node_id++);
511  }
512  }
513  else if (ix == _nx - 1 && iy == _ny - 1) // top right corner
514  {
515  for (auto z : _z_grid)
516  {
517  Point z0{0, 0, z};
518  for (unsigned int i = 0; i < points_per_corner; i++)
519  mesh_base->add_point(tr_corner_points[i] + x0 + y0 + z0, node_id++);
520  }
521  }
522  else if (ix == 0 && (iy != _ny - 1 || iy != 0)) // left side
523  {
524  for (auto z : _z_grid)
525  {
526  Point z0{0, 0, z};
527  for (unsigned int i = 0; i < points_per_side; i++)
528  mesh_base->add_point(left_points[i] + x0 + y0 + z0, node_id++);
529  }
530  }
531  else if (ix == _nx - 1 && (iy != _ny - 1 || iy != 0)) // right side
532  {
533  for (auto z : _z_grid)
534  {
535  Point z0{0, 0, z};
536  for (unsigned int i = 0; i < points_per_side; i++)
537  mesh_base->add_point(right_points[i] + x0 + y0 + z0, node_id++);
538  }
539  }
540  else if (iy == 0 && (ix != _nx - 1 || ix != 0)) // bottom side
541  {
542  for (auto z : _z_grid)
543  {
544  Point z0{0, 0, z};
545  for (unsigned int i = 0; i < points_per_side; i++)
546  mesh_base->add_point(bottom_points[i] + x0 + y0 + z0, node_id++);
547  }
548  }
549  else if (iy == _ny - 1 && (ix != _nx - 1 || ix != 0)) // top side
550  {
551  for (auto z : _z_grid)
552  {
553  Point z0{0, 0, z};
554  for (unsigned int i = 0; i < points_per_side; i++)
555  mesh_base->add_point(top_points[i] + x0 + y0 + z0, node_id++);
556  }
557  }
558  else // center
559  {
560  for (auto z : _z_grid)
561  {
562  Point z0{0, 0, z};
563  for (unsigned int i = 0; i < points_per_center; i++)
564  mesh_base->add_point(center_points[i] + x0 + y0 + z0, node_id++);
565  }
566  }
567  }
568  }
569  }
570 
571  // Add elements to the mesh. The elements are 6-node prisms. The
572  // bases of these prisms form a triangulated representation of a cross-section
573  // of a center subchannel.
574  if (_n_channels == 2)
575  {
576  unsigned int elem_id = 0;
577  for (unsigned int iy = 0; iy < _ny; iy++)
578  {
579  for (unsigned int ix = 0; ix < _nx; ix++)
580  {
581  unsigned int i_ch = _nx * iy + ix;
582  for (unsigned int iz = 0; iz < _n_cells; iz++)
583  {
584  for (unsigned int i = 0; i < elems_per_side; i++)
585  {
586  Elem * elem = new Prism6;
587  elem->subdomain_id() = _block_id;
588  elem->set_id(elem_id++);
589  elem = mesh_base->add_elem(elem);
590  // index of the central node at base of cell
591  unsigned int indx1 = iz * points_per_side + points_per_side * (_n_cells + 1) * i_ch;
592  // index of the central node at top of cell
593  unsigned int indx2 =
594  (iz + 1) * points_per_side + points_per_side * (_n_cells + 1) * i_ch;
595  unsigned int elems_per_channel = elems_per_side;
596  elem->set_node(0, mesh_base->node_ptr(indx1));
597  elem->set_node(1, mesh_base->node_ptr(indx1 + i + 1));
598  if (i != elems_per_channel - 1)
599  elem->set_node(2, mesh_base->node_ptr(indx1 + i + 2));
600  else
601  elem->set_node(2, mesh_base->node_ptr(indx1 + 1));
602 
603  elem->set_node(3, mesh_base->node_ptr(indx2));
604  elem->set_node(4, mesh_base->node_ptr(indx2 + i + 1));
605  if (i != elems_per_channel - 1)
606  elem->set_node(5, mesh_base->node_ptr(indx2 + i + 2));
607  else
608  elem->set_node(5, mesh_base->node_ptr(indx2 + 1));
609 
610  if (iz == 0)
611  boundary_info.add_side(elem, 0, 0);
612  if (iz == _n_cells - 1)
613  boundary_info.add_side(elem, 4, 1);
614  }
615  }
616  }
617  }
618  boundary_info.sideset_name(0) = "inlet";
619  boundary_info.sideset_name(1) = "outlet";
620  mesh_base->subdomain_name(_block_id) = name();
621  mesh_base->prepare_for_use();
622  }
623  else if (_n_channels > 2 && (_ny == 1 || _nx == 1))
624  {
625  unsigned int elem_id = 0;
626  unsigned int number_of_corner = 0;
627  unsigned int number_of_side = 0;
628  unsigned int number_of_center = 0;
629  unsigned int elems_per_channel = 0;
630  unsigned int points_per_channel = 0;
631  for (unsigned int iy = 0; iy < _ny; iy++)
632  {
633  for (unsigned int ix = 0; ix < _nx; ix++)
634  {
635  unsigned int i_ch = _nx * iy + ix;
636  auto subch_type = getSubchannelType(i_ch);
637  // note that in this case i use side geometry for corner subchannel
638  if (subch_type == EChannelType::CORNER)
639  {
640  number_of_side++;
641  elems_per_channel = elems_per_side;
642  points_per_channel = points_per_side;
643  }
644  // note that in this case i use center geometry for edge subchannel
645  else if (subch_type == EChannelType::EDGE)
646  {
647  number_of_center++;
648  elems_per_channel = elems_per_center;
649  points_per_channel = points_per_center;
650  }
651  for (unsigned int iz = 0; iz < _n_cells; iz++)
652  {
653  unsigned int elapsed_points = number_of_corner * points_per_corner * (_n_cells + 1) +
654  number_of_side * points_per_side * (_n_cells + 1) +
655  number_of_center * points_per_center * (_n_cells + 1) -
656  points_per_channel * (_n_cells + 1);
657  // index of the central node at base of cell
658  unsigned int indx1 = iz * points_per_channel + elapsed_points;
659  // index of the central node at top of cell
660  unsigned int indx2 = (iz + 1) * points_per_channel + elapsed_points;
661 
662  for (unsigned int i = 0; i < elems_per_channel; i++)
663  {
664  Elem * elem = new Prism6;
665  elem->subdomain_id() = _block_id;
666  elem->set_id(elem_id++);
667  elem = mesh_base->add_elem(elem);
668 
669  elem->set_node(0, mesh_base->node_ptr(indx1));
670  elem->set_node(1, mesh_base->node_ptr(indx1 + i + 1));
671  if (i != elems_per_channel - 1)
672  elem->set_node(2, mesh_base->node_ptr(indx1 + i + 2));
673  else
674  elem->set_node(2, mesh_base->node_ptr(indx1 + 1));
675 
676  elem->set_node(3, mesh_base->node_ptr(indx2));
677  elem->set_node(4, mesh_base->node_ptr(indx2 + i + 1));
678  if (i != elems_per_channel - 1)
679  elem->set_node(5, mesh_base->node_ptr(indx2 + i + 2));
680  else
681  elem->set_node(5, mesh_base->node_ptr(indx2 + 1));
682 
683  if (iz == 0)
684  boundary_info.add_side(elem, 0, 0);
685  if (iz == _n_cells - 1)
686  boundary_info.add_side(elem, 4, 1);
687  }
688  }
689  }
690  }
691  boundary_info.sideset_name(0) = "inlet";
692  boundary_info.sideset_name(1) = "outlet";
693  mesh_base->subdomain_name(_block_id) = name();
694  mesh_base->prepare_for_use();
695  }
696  else
697  {
698  unsigned int elem_id = 0;
699  unsigned int number_of_corner = 0;
700  unsigned int number_of_side = 0;
701  unsigned int number_of_center = 0;
702  unsigned int elems_per_channel = 0;
703  unsigned int points_per_channel = 0;
704  for (unsigned int iy = 0; iy < _ny; iy++)
705  {
706  for (unsigned int ix = 0; ix < _nx; ix++)
707  {
708  unsigned int i_ch = _nx * iy + ix;
709  auto subch_type = getSubchannelType(i_ch);
710  if (subch_type == EChannelType::CORNER)
711  {
712  number_of_corner++;
713  elems_per_channel = elems_per_corner;
714  points_per_channel = points_per_corner;
715  }
716  else if (subch_type == EChannelType::EDGE)
717  {
718  number_of_side++;
719  elems_per_channel = elems_per_side;
720  points_per_channel = points_per_side;
721  }
722  else
723  {
724  number_of_center++;
725  elems_per_channel = elems_per_center;
726  points_per_channel = points_per_center;
727  }
728  for (unsigned int iz = 0; iz < _n_cells; iz++)
729  {
730  unsigned int elapsed_points = number_of_corner * points_per_corner * (_n_cells + 1) +
731  number_of_side * points_per_side * (_n_cells + 1) +
732  number_of_center * points_per_center * (_n_cells + 1) -
733  points_per_channel * (_n_cells + 1);
734  // index of the central node at base of cell
735  unsigned int indx1 = iz * points_per_channel + elapsed_points;
736  // index of the central node at top of cell
737  unsigned int indx2 = (iz + 1) * points_per_channel + elapsed_points;
738 
739  for (unsigned int i = 0; i < elems_per_channel; i++)
740  {
741  Elem * elem = new Prism6;
742  elem->subdomain_id() = _block_id;
743  elem->set_id(elem_id++);
744  elem = mesh_base->add_elem(elem);
745 
746  elem->set_node(0, mesh_base->node_ptr(indx1));
747  elem->set_node(1, mesh_base->node_ptr(indx1 + i + 1));
748  if (i != elems_per_channel - 1)
749  elem->set_node(2, mesh_base->node_ptr(indx1 + i + 2));
750  else
751  elem->set_node(2, mesh_base->node_ptr(indx1 + 1));
752 
753  elem->set_node(3, mesh_base->node_ptr(indx2));
754  elem->set_node(4, mesh_base->node_ptr(indx2 + i + 1));
755  if (i != elems_per_channel - 1)
756  elem->set_node(5, mesh_base->node_ptr(indx2 + i + 2));
757  else
758  elem->set_node(5, mesh_base->node_ptr(indx2 + 1));
759 
760  if (iz == 0)
761  boundary_info.add_side(elem, 0, 0);
762  if (iz == _n_cells - 1)
763  boundary_info.add_side(elem, 4, 1);
764  }
765  }
766  }
767  }
768  boundary_info.sideset_name(0) = "inlet";
769  boundary_info.sideset_name(1) = "outlet";
770  mesh_base->subdomain_name(_block_id) = name();
771  mesh_base->prepare_for_use();
772  }
773 
774  return mesh_base;
775 }
const Real radius
const Real _pitch
Distance between the neighbor fuel pins, pitch.
const unsigned int _n_cells
Number of cells in the axial direction.
std::vector< Real > _z_grid
axial location of nodes
unsigned int _n_channels
Total number of subchannels.
const unsigned int _nx
Number of subchannels in the x direction.
const std::string & name() const
const Real _side_gap
The side gap, not to be confused with the gap between pins, this refers to the gap next to the duct o...
const unsigned int _ny
Number of subchannels in the y direction.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
EChannelType getSubchannelType(unsigned int index) const
std::unique_ptr< MeshBase > buildMeshBaseObject(unsigned int dim=libMesh::invalid_uint)
const unsigned int & _block_id
Subdomain ID used for the mesh block.

◆ getSubchannelPosition()

std::vector<Real> SCMDetailedQuadSubChannelMeshGenerator::getSubchannelPosition ( unsigned int  i)
inlineprotected

Definition at line 26 of file SCMDetailedQuadSubChannelMeshGenerator.h.

26 { return _subchannel_position[i]; }
std::vector< std::vector< Real > > _subchannel_position
x,y coordinates of the subchannel centroids

◆ getSubchannelType()

EChannelType SCMDetailedQuadSubChannelMeshGenerator::getSubchannelType ( unsigned int  index) const
inlineprotected

Definition at line 25 of file SCMDetailedQuadSubChannelMeshGenerator.h.

Referenced by generate().

25 { return _subch_type[index]; }
std::vector< EChannelType > _subch_type
Subchannel type.

◆ validParams()

InputParameters SCMDetailedQuadSubChannelMeshGenerator::validParams ( )
static

Definition at line 23 of file SCMDetailedQuadSubChannelMeshGenerator.C.

24 {
26  params.addClassDescription(
27  "Creates a detailed mesh of subchannels in a square lattice arrangement");
28  params.addRequiredParam<Real>("pitch", "Pitch [m]");
29  params.addRequiredParam<Real>("pin_diameter", "Rod diameter [m]");
30  params.addParam<Real>("unheated_length_entry", 0.0, "Unheated length at entry [m]");
31  params.addRequiredParam<Real>("heated_length", "Heated length [m]");
32  params.addParam<Real>("unheated_length_exit", 0.0, "Unheated length at exit [m]");
33  params.addRequiredParam<unsigned int>("n_cells", "The number of cells in the axial direction");
34  params.addRequiredParam<unsigned int>("nx", "Number of channels in the x direction [-]");
35  params.addRequiredParam<unsigned int>("ny", "Number of channels in the y direction [-]");
36  params.addRequiredParam<Real>(
37  "gap",
38  "The side gap, not to be confused with the gap between pins, this refers to the gap "
39  "next to the duct or else the distance between the subchannel centroid to the duct wall."
40  "distance(edge pin center, duct wall) = pitch / 2 + side_gap [m]");
41  params.deprecateParam("gap", "side_gap", "08/06/2026");
42  params.addParam<unsigned int>("block_id", 0, "Block ID used for the mesh subdomain.");
43  return params;
44 }
void addParam(const std::string &name, const std::initializer_list< typename T::value_type > &value, const std::string &doc_string)
void addRequiredParam(const std::string &name, const std::string &doc_string)
void deprecateParam(const std::string &old_name, const std::string &new_name, const std::string &removal_date)
static InputParameters validParams()
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void addClassDescription(const std::string &doc_string)

Member Data Documentation

◆ _block_id

const unsigned int& SCMDetailedQuadSubChannelMeshGenerator::_block_id
protected

Subdomain ID used for the mesh block.

Definition at line 59 of file SCMDetailedQuadSubChannelMeshGenerator.h.

Referenced by generate().

◆ _heated_length

const Real SCMDetailedQuadSubChannelMeshGenerator::_heated_length
protected

heated length of the fuel Pin

Definition at line 31 of file SCMDetailedQuadSubChannelMeshGenerator.h.

Referenced by SCMDetailedQuadSubChannelMeshGenerator().

◆ _n_cells

const unsigned int SCMDetailedQuadSubChannelMeshGenerator::_n_cells
protected

Number of cells in the axial direction.

Definition at line 41 of file SCMDetailedQuadSubChannelMeshGenerator.h.

Referenced by generate(), and SCMDetailedQuadSubChannelMeshGenerator().

◆ _n_channels

unsigned int SCMDetailedQuadSubChannelMeshGenerator::_n_channels
protected

Total number of subchannels.

Definition at line 47 of file SCMDetailedQuadSubChannelMeshGenerator.h.

Referenced by generate(), and SCMDetailedQuadSubChannelMeshGenerator().

◆ _nx

const unsigned int SCMDetailedQuadSubChannelMeshGenerator::_nx
protected

Number of subchannels in the x direction.

Definition at line 43 of file SCMDetailedQuadSubChannelMeshGenerator.h.

Referenced by generate(), and SCMDetailedQuadSubChannelMeshGenerator().

◆ _ny

const unsigned int SCMDetailedQuadSubChannelMeshGenerator::_ny
protected

Number of subchannels in the y direction.

Definition at line 45 of file SCMDetailedQuadSubChannelMeshGenerator.h.

Referenced by generate(), and SCMDetailedQuadSubChannelMeshGenerator().

◆ _pin_diameter

const Real SCMDetailedQuadSubChannelMeshGenerator::_pin_diameter
protected

fuel Pin diameter

Definition at line 39 of file SCMDetailedQuadSubChannelMeshGenerator.h.

Referenced by generate().

◆ _pitch

const Real SCMDetailedQuadSubChannelMeshGenerator::_pitch
protected

Distance between the neighbor fuel pins, pitch.

Definition at line 37 of file SCMDetailedQuadSubChannelMeshGenerator.h.

Referenced by generate(), and SCMDetailedQuadSubChannelMeshGenerator().

◆ _side_gap

const Real SCMDetailedQuadSubChannelMeshGenerator::_side_gap
protected

The side gap, not to be confused with the gap between pins, this refers to the gap next to the duct or else the distance between the subchannel centroid to the duct wall.

distance(edge pin center, duct wall) = pitch / 2 + side_gap [m].

Definition at line 53 of file SCMDetailedQuadSubChannelMeshGenerator.h.

Referenced by generate().

◆ _subch_type

std::vector<EChannelType> SCMDetailedQuadSubChannelMeshGenerator::_subch_type
protected

Subchannel type.

Definition at line 55 of file SCMDetailedQuadSubChannelMeshGenerator.h.

Referenced by getSubchannelType(), and SCMDetailedQuadSubChannelMeshGenerator().

◆ _subchannel_position

std::vector<std::vector<Real> > SCMDetailedQuadSubChannelMeshGenerator::_subchannel_position
protected

x,y coordinates of the subchannel centroids

Definition at line 57 of file SCMDetailedQuadSubChannelMeshGenerator.h.

Referenced by getSubchannelPosition(), and SCMDetailedQuadSubChannelMeshGenerator().

◆ _unheated_length_entry

const Real SCMDetailedQuadSubChannelMeshGenerator::_unheated_length_entry
protected

unheated length of the fuel Pin at the entry of the assembly

Definition at line 29 of file SCMDetailedQuadSubChannelMeshGenerator.h.

Referenced by SCMDetailedQuadSubChannelMeshGenerator().

◆ _unheated_length_exit

const Real SCMDetailedQuadSubChannelMeshGenerator::_unheated_length_exit
protected

unheated length of the fuel Pin at the exit of the assembly

Definition at line 33 of file SCMDetailedQuadSubChannelMeshGenerator.h.

Referenced by SCMDetailedQuadSubChannelMeshGenerator().

◆ _z_grid

std::vector<Real> SCMDetailedQuadSubChannelMeshGenerator::_z_grid
protected

axial location of nodes

Definition at line 35 of file SCMDetailedQuadSubChannelMeshGenerator.h.

Referenced by generate(), and SCMDetailedQuadSubChannelMeshGenerator().


The documentation for this class was generated from the following files: