Loading [MathJax]/extensions/tex2jax.js
https://mooseframework.inl.gov
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends
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
 
virtual const std::string & name () const
 
std::string typeAndName () const
 
std::string errorPrefix (const std::string &error_type) const
 
void callMooseError (std::string msg, const bool with_prefix) const
 
MooseObjectParameterName uniqueParameterName (const std::string &parameter_name) const
 
const InputParametersparameters () const
 
MooseObjectName uniqueName () 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 & 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 &nm) 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
 
void connectControllableParams (const std::string &parameter, const std::string &object_type, const std::string &object_name, const std::string &object_parameter) const
 
void mooseError (Args &&... args) const
 
void mooseErrorNonPrefixed (Args &&... args) const
 
void mooseDocumentedError (const std::string &repo_name, const unsigned int issue_num, Args &&... args) const
 
void mooseWarning (Args &&... args) const
 
void mooseWarningNonPrefixed (Args &&... args) const
 
void mooseDeprecated (Args &&... args) const
 
void mooseInfo (Args &&... args) 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 setHasGenerateData (InputParameters &params)
 

Public Attributes

const ConsoleStream _console
 

Static Public Attributes

static const std::string data_only_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 _gap
 Half of gap between adjacent assemblies. 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
 
const std::string _type
 
const std::string _name
 
const InputParameters_pars
 
Factory_factory
 
ActionFactory_action_factory
 
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 43 of file SCMDetailedQuadSubChannelMeshGenerator.C.

46  _unheated_length_entry(getParam<Real>("unheated_length_entry")),
47  _heated_length(getParam<Real>("heated_length")),
48  _unheated_length_exit(getParam<Real>("unheated_length_exit")),
49  _pitch(getParam<Real>("pitch")),
50  _pin_diameter(getParam<Real>("pin_diameter")),
51  _n_cells(getParam<unsigned int>("n_cells")),
52  _nx(getParam<unsigned int>("nx")),
53  _ny(getParam<unsigned int>("ny")),
54  _n_channels(_nx * _ny),
55  _gap(getParam<Real>("gap")),
56  _block_id(getParam<unsigned int>("block_id"))
57 {
59  Real dz = L / _n_cells;
60  for (unsigned int i = 0; i < _n_cells + 1; i++)
61  _z_grid.push_back(dz * i);
62 
64  for (unsigned int i = 0; i < _n_channels; i++)
65  {
66  _subchannel_position[i].reserve(3);
67  for (unsigned int j = 0; j < 3; j++)
68  {
69  _subchannel_position.at(i).push_back(0.0);
70  }
71  }
72 
73  _subch_type.resize(_n_channels);
74  for (unsigned int iy = 0; iy < _ny; iy++)
75  {
76  for (unsigned int ix = 0; ix < _nx; ix++)
77  {
78  unsigned int i_ch = _nx * iy + ix;
79  bool is_corner = (ix == 0 && iy == 0) || (ix == _nx - 1 && iy == 0) ||
80  (ix == 0 && iy == _ny - 1) || (ix == _nx - 1 && iy == _ny - 1);
81  bool is_edge = (ix == 0 || iy == 0 || ix == _nx - 1 || iy == _ny - 1);
82 
83  if (is_corner)
85  else if (is_edge)
87  else
89 
90  // set the subchannel positions so that the center of the assembly is the zero point
91  Real offset_x = (_nx - 1) * _pitch / 2.0;
92  Real offset_y = (_ny - 1) * _pitch / 2.0;
93  _subchannel_position[i_ch][0] = _pitch * ix - offset_x;
94  _subchannel_position[i_ch][1] = _pitch * iy - offset_y;
95  }
96  }
97 }
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 Real _gap
Half of gap between adjacent assemblies.
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 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.
const InputParameters & parameters() const
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 100 of file SCMDetailedQuadSubChannelMeshGenerator.C.

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

Referenced by generate().

◆ _gap

const Real SCMDetailedQuadSubChannelMeshGenerator::_gap
protected

Half of gap between adjacent assemblies.

Definition at line 49 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().

◆ _subch_type

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

Subchannel type.

Definition at line 51 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 53 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: