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 "SCMTriAssemblyMeshGenerator.h"
11 : #include "TriSubChannelMesh.h"
12 : #include <cmath>
13 : #include <memory>
14 : #include "libmesh/edge_edge2.h"
15 : #include "libmesh/unstructured_mesh.h"
16 :
17 : registerMooseObject("SubChannelApp", SCMTriAssemblyMeshGenerator);
18 : registerMooseObjectRenamed("SubChannelApp",
19 : SCMTriSubChannelMeshGenerator,
20 : "06/30/2027 24:00",
21 : SCMTriAssemblyMeshGenerator);
22 : registerMooseObjectRenamed("SubChannelApp",
23 : TriSubChannelMeshGenerator,
24 : "06/30/2027 24:00",
25 : SCMTriAssemblyMeshGenerator);
26 : registerMooseObjectRenamed("SubChannelApp",
27 : SCMTriPinMeshGenerator,
28 : "06/30/2027 24:00",
29 : SCMTriAssemblyMeshGenerator);
30 : registerMooseObjectRenamed("SubChannelApp",
31 : TriPinMeshGenerator,
32 : "06/30/2027 24:00",
33 : SCMTriAssemblyMeshGenerator);
34 :
35 : InputParameters
36 368 : SCMTriAssemblyMeshGenerator::validParams()
37 : {
38 368 : InputParameters params = MeshGenerator::validParams();
39 368 : params.addClassDescription(
40 : "Creates a mesh of 1D subchannels and 1D pins in a triangular lattice arrangement");
41 736 : params.addRequiredParam<unsigned int>("n_cells", "The number of cells in the axial direction");
42 736 : params.addRequiredParam<Real>("pitch", "Pitch [m]");
43 736 : params.addRequiredParam<Real>("pin_diameter", "Rod diameter [m]");
44 736 : params.addParam<Real>("unheated_length_entry", 0.0, "Unheated length at entry [m]");
45 736 : params.addRequiredParam<Real>("heated_length", "Heated length [m]");
46 736 : params.addParam<Real>("unheated_length_exit", 0.0, "Unheated length at exit [m]");
47 736 : params.addRequiredParam<unsigned int>(
48 : "nrings",
49 : "Number of fuel-pin rings per assembly, counting the center pin as the first ring [-]");
50 736 : params.addRequiredParam<Real>("flat_to_flat",
51 : "Flat to flat distance for the hexagonal assembly [m]");
52 736 : params.addRequiredParam<Real>("dwire", "Wire diameter [m]");
53 736 : params.addRequiredParam<Real>("hwire", "Wire lead length [m]");
54 736 : params.addParam<std::vector<Real>>(
55 : "spacer_z", {}, "Axial location of spacers/vanes/mixing vanes [m]");
56 736 : params.addParam<std::vector<Real>>(
57 : "spacer_k", {}, "K-loss coefficient of spacers/vanes/mixing vanes [-]");
58 736 : params.addParam<Real>("Kij", 0.5, "Lateral form loss coefficient [-]");
59 736 : params.addParam<std::vector<Real>>("z_blockage",
60 736 : std::vector<Real>({0.0, 0.0}),
61 : "axial location of blockage (inlet, outlet) [m]");
62 736 : params.addParam<std::vector<unsigned int>>("index_blockage",
63 736 : std::vector<unsigned int>({0}),
64 : "index of subchannels affected by blockage");
65 736 : params.addParam<std::vector<Real>>(
66 : "reduction_blockage",
67 736 : std::vector<Real>({1.0}),
68 : "Area reduction of subchannels affected by blockage (number to muliply the area)");
69 736 : params.addParam<std::vector<Real>>("k_blockage",
70 736 : std::vector<Real>({0.0}),
71 : "Form loss coefficient of subchannels affected by blockage");
72 736 : params.addParam<unsigned int>("block_id", 0, "Subchannel block id");
73 736 : params.deprecateParam("block_id", "subchannel_block_id", "07/01/2027");
74 736 : params.addParam<unsigned int>("pin_block_id", 1, "Fuel Pin block id");
75 368 : return params;
76 0 : }
77 :
78 184 : SCMTriAssemblyMeshGenerator::SCMTriAssemblyMeshGenerator(const InputParameters & params)
79 : : MeshGenerator(params),
80 184 : _unheated_length_entry(getParam<Real>("unheated_length_entry")),
81 368 : _heated_length(getParam<Real>("heated_length")),
82 368 : _unheated_length_exit(getParam<Real>("unheated_length_exit")),
83 368 : _subchannel_block_id(getParam<unsigned int>("subchannel_block_id")),
84 368 : _pin_block_id(getParam<unsigned int>("pin_block_id")),
85 368 : _spacer_z(getParam<std::vector<Real>>("spacer_z")),
86 368 : _spacer_k(getParam<std::vector<Real>>("spacer_k")),
87 368 : _z_blockage(getParam<std::vector<Real>>("z_blockage")),
88 368 : _index_blockage(getParam<std::vector<unsigned int>>("index_blockage")),
89 368 : _reduction_blockage(getParam<std::vector<Real>>("reduction_blockage")),
90 368 : _k_blockage(getParam<std::vector<Real>>("k_blockage")),
91 368 : _pitch(getParam<Real>("pitch")),
92 368 : _kij(getParam<Real>("Kij")),
93 368 : _pin_diameter(getParam<Real>("pin_diameter")),
94 368 : _n_cells(getParam<unsigned int>("n_cells")),
95 368 : _n_rings(getParam<unsigned int>("nrings")),
96 184 : _n_channels(0),
97 368 : _flat_to_flat(getParam<Real>("flat_to_flat")),
98 368 : _dwire(getParam<Real>("dwire")),
99 368 : _hwire(getParam<Real>("hwire")),
100 184 : _duct_to_pin_gap(0.5 *
101 184 : (_flat_to_flat - (_n_rings - 1) * _pitch * std::sqrt(3.0) - _pin_diameter)),
102 184 : _npins(0),
103 184 : _n_gaps(0)
104 : {
105 184 : const Real total_length = _unheated_length_entry + _heated_length + _unheated_length_exit;
106 :
107 184 : if (_n_rings < 2)
108 0 : paramError("nrings",
109 : "'nrings' must be at least 2. In this mesh generator, the center pin counts as "
110 : "the first ring, so a 7-pin bundle uses nrings = 2.");
111 :
112 184 : if (_n_cells == 0)
113 0 : paramError("n_cells", "The number of axial cells must be greater than zero");
114 :
115 184 : if (total_length <= 0.0)
116 0 : mooseError("Total bundle length must be greater than zero");
117 :
118 184 : if (_spacer_z.size() != _spacer_k.size())
119 0 : mooseError("Size of vector spacer_z should be equal to size of vector spacer_k");
120 :
121 354 : for (const auto spacer_z : _spacer_z)
122 170 : if (spacer_z < 0.0 || spacer_z > total_length)
123 0 : paramError("spacer_z", "Location of spacers should be between zero and total bundle length");
124 :
125 184 : if (_z_blockage.size() != 2)
126 0 : paramError("z_blockage", "Size of vector z_blockage must be 2");
127 :
128 184 : if (_z_blockage.front() > _z_blockage.back())
129 0 : paramError("z_blockage", "z_blockage inlet location must not exceed outlet location");
130 :
131 184 : if (*max_element(_reduction_blockage.begin(), _reduction_blockage.end()) > 1)
132 0 : paramError("reduction_blockage",
133 : "The area reduction of the blocked subchannels cannot be more than 1");
134 :
135 184 : if ((_index_blockage.size() != _reduction_blockage.size()) ||
136 368 : (_index_blockage.size() != _k_blockage.size()) ||
137 : (_reduction_blockage.size() != _k_blockage.size()))
138 0 : mooseError("Size of vectors: index_blockage, reduction_blockage, k_blockage, must be equal "
139 : "to eachother");
140 :
141 184 : SubChannelMesh::generateZGrid(
142 184 : _unheated_length_entry, _heated_length, _unheated_length_exit, _n_cells, _z_grid);
143 :
144 : // Defining the total length from 3 axial sections
145 : Real L = total_length;
146 :
147 : // Defining the position of the spacer grid in the numerical solution array
148 : std::vector<int> spacer_cell;
149 354 : for (const auto & elem : _spacer_z)
150 170 : spacer_cell.emplace_back(std::round(elem * _n_cells / L));
151 :
152 : // Defining the array for axial resistances
153 : std::vector<Real> kgrid;
154 184 : kgrid.resize(_n_cells + 1, 0.0);
155 :
156 : // Summing the spacer resistance to the grid resistance array
157 354 : for (unsigned int index = 0; index < spacer_cell.size(); index++)
158 170 : kgrid[spacer_cell[index]] += _spacer_k[index];
159 :
160 : // compute the hex mesh variables
161 : // -------------------------------------------
162 :
163 : // x coordinate for the first position
164 : Real x0 = 0.0;
165 : // y coordinate for the first position
166 : Real y0 = 0.0;
167 : // x coordinate for the second position
168 : Real x1 = 0.0;
169 : // y coordinate for the second position dummy variable
170 : Real y1 = 0.0;
171 : // dummy variable
172 : Real a1 = 0.0;
173 : // dummy variable
174 : Real a2 = 0.0;
175 : // average x coordinate
176 : Real avg_coor_x = 0.0;
177 : // average y coordinate
178 : Real avg_coor_y = 0.0;
179 : // distance between two points
180 : Real dist = 0.0;
181 : // distance between two points
182 : Real dist0 = 0.0;
183 : // integer counter
184 184 : unsigned int kgap = 0;
185 : // dummy integer
186 : unsigned int icorner = 0;
187 : // used to defined global direction of the cross_flow_map coefficients for each subchannel and gap
188 184 : const Real positive_flow = 1.0;
189 : // used to defined global direction of the cross_flow_map coefficients for each subchannel and gap
190 : const Real negative_flow = -1.0;
191 : // the indicator used while setting _gap_to_chan_map array
192 : std::vector<std::pair<unsigned int, unsigned int>> gap_fill;
193 184 : TriSubChannelMesh::pinPositions(_pin_position, _n_rings, _pitch, Point(0, 0));
194 184 : _npins = _pin_position.size();
195 : // assign the pins to the corresponding rings
196 : unsigned int k = 0; // initialize the fuel Pin counter index
197 184 : _pins_in_rings.resize(_n_rings);
198 184 : _pins_in_rings[0].push_back(k++);
199 727 : for (unsigned int i = 1; i < _n_rings; i++)
200 8097 : for (unsigned int j = 0; j < i * 6; j++)
201 7554 : _pins_in_rings[i].push_back(k++);
202 : // Given the number of pins and number of fuel Pin rings, the number of subchannels can be
203 : // computed as follows:
204 : unsigned int chancount = 0.0;
205 : // Summing internal channels
206 727 : for (unsigned int j = 0; j < _n_rings - 1; j++)
207 543 : chancount += j * 6;
208 : // Adding external channels to the total count
209 184 : _n_channels = chancount + _npins - 1 + (_n_rings - 1) * 6 + 6;
210 :
211 184 : if (*max_element(_index_blockage.begin(), _index_blockage.end()) > (_n_channels - 1))
212 0 : paramError("index_blockage",
213 : "The index of the blocked subchannel cannot be more than the max index of the "
214 : "subchannels");
215 :
216 184 : if ((_index_blockage.size() > _n_channels) || (_reduction_blockage.size() > _n_channels) ||
217 : (_k_blockage.size() > _n_channels))
218 0 : mooseError("Size of vectors: index_blockage, reduction_blockage, k_blockage, cannot be more "
219 : "than the total number of subchannels");
220 :
221 : // Defining the 2D array for axial resistances
222 184 : _k_grid.resize(_n_channels, std::vector<Real>(_n_cells + 1));
223 16396 : for (unsigned int i = 0; i < _n_channels; i++)
224 16212 : _k_grid[i] = kgrid;
225 :
226 : // Add blockage resistance to the 2D grid resistane array
227 184 : Real dz = L / _n_cells;
228 4705 : for (unsigned int i = 0; i < _n_cells + 1; i++)
229 : {
230 4521 : if ((dz * i >= _z_blockage.front() && dz * i <= _z_blockage.back()))
231 : {
232 : unsigned int index(0);
233 585 : for (const auto & i_ch : _index_blockage)
234 : {
235 401 : _k_grid[i_ch][i] += _k_blockage[index];
236 401 : index++;
237 : }
238 : }
239 : }
240 :
241 184 : _chan_to_pin_map.resize(_n_channels);
242 184 : _pin_to_chan_map.resize(_npins);
243 184 : _subch_type.resize(_n_channels);
244 184 : _n_gaps = _n_channels + _npins - 1; /// initial assignment
245 184 : _gap_to_chan_map.resize(_n_gaps);
246 184 : _gap_to_pin_map.resize(_n_gaps);
247 184 : gap_fill.resize(_n_gaps);
248 184 : _chan_to_gap_map.resize(_n_channels);
249 184 : _gap_pairs_sf.resize(_n_channels);
250 184 : _chan_pairs_sf.resize(_n_channels);
251 184 : _gij_map.resize(_n_cells + 1);
252 184 : _sign_id_crossflow_map.resize(_n_channels);
253 184 : _gap_type.resize(_n_gaps);
254 184 : _subchannel_position.resize(_n_channels);
255 :
256 16396 : for (unsigned int i = 0; i < _n_channels; i++)
257 : {
258 16212 : _chan_to_pin_map[i].reserve(3);
259 16212 : _chan_to_gap_map[i].reserve(3);
260 16212 : _sign_id_crossflow_map[i].reserve(3);
261 16212 : _subchannel_position[i].reserve(3);
262 64848 : for (unsigned int j = 0; j < 3; j++)
263 : {
264 48636 : _sign_id_crossflow_map.at(i).push_back(positive_flow);
265 48636 : _subchannel_position.at(i).push_back(0.0);
266 : }
267 : }
268 :
269 4705 : for (unsigned int iz = 0; iz < _n_cells + 1; iz++)
270 : {
271 4521 : _gij_map[iz].reserve(_n_gaps);
272 : }
273 :
274 7922 : for (unsigned int i = 0; i < _npins; i++)
275 7738 : _pin_to_chan_map[i].reserve(6);
276 :
277 : // create the subchannels
278 : k = 0; // initialize the subchannel counter index
279 : kgap = 0;
280 : // for each ring we trace the subchannels by pairing up to neighbor pins and looking for the third
281 : // Pin at inner or outer ring compared to the current ring.
282 727 : for (unsigned int i = 1; i < _n_rings; i++)
283 : {
284 : // find the closest Pin at back ring
285 8097 : for (unsigned int j = 0; j < _pins_in_rings[i].size(); j++)
286 : {
287 7554 : if (j == _pins_in_rings[i].size() - 1)
288 : {
289 543 : _chan_to_pin_map[k].push_back(_pins_in_rings[i][j]);
290 543 : _chan_to_pin_map[k].push_back(_pins_in_rings[i][0]);
291 543 : avg_coor_x =
292 543 : 0.5 * (_pin_position[_pins_in_rings[i][j]](0) + _pin_position[_pins_in_rings[i][0]](0));
293 543 : avg_coor_y =
294 543 : 0.5 * (_pin_position[_pins_in_rings[i][j]](1) + _pin_position[_pins_in_rings[i][0]](1));
295 543 : _gap_to_pin_map[kgap].first = _pins_in_rings[i][0];
296 543 : _gap_to_pin_map[kgap].second = _pins_in_rings[i][j];
297 543 : _gap_type[kgap] = EChannelType::CENTER;
298 543 : kgap = kgap + 1;
299 : }
300 : else
301 : {
302 7011 : _chan_to_pin_map[k].push_back(_pins_in_rings[i][j]);
303 7011 : _chan_to_pin_map[k].push_back(_pins_in_rings[i][j + 1]);
304 7011 : avg_coor_x = 0.5 * (_pin_position[_pins_in_rings[i][j]](0) +
305 7011 : _pin_position[_pins_in_rings[i][j + 1]](0));
306 7011 : avg_coor_y = 0.5 * (_pin_position[_pins_in_rings[i][j]](1) +
307 7011 : _pin_position[_pins_in_rings[i][j + 1]](1));
308 7011 : _gap_to_pin_map[kgap].first = _pins_in_rings[i][j];
309 7011 : _gap_to_pin_map[kgap].second = _pins_in_rings[i][j + 1];
310 7011 : _gap_type[kgap] = EChannelType::CENTER;
311 7011 : kgap = kgap + 1;
312 : }
313 :
314 : dist0 = 1.0e+5;
315 :
316 7554 : _chan_to_pin_map[k].push_back(_pins_in_rings[i - 1][0]);
317 : unsigned int l0 = 0;
318 :
319 109818 : for (unsigned int l = 0; l < _pins_in_rings[i - 1].size(); l++)
320 : {
321 102264 : dist = std::sqrt(pow(_pin_position[_pins_in_rings[i - 1][l]](0) - avg_coor_x, 2) +
322 102264 : pow(_pin_position[_pins_in_rings[i - 1][l]](1) - avg_coor_y, 2));
323 :
324 102264 : if (dist < dist0)
325 : {
326 34701 : _chan_to_pin_map[k][2] = _pins_in_rings[i - 1][l];
327 : l0 = l;
328 : dist0 = dist;
329 : } // if
330 : } // l
331 :
332 7554 : _gap_to_pin_map[kgap].first = _pins_in_rings[i][j];
333 7554 : _gap_to_pin_map[kgap].second = _pins_in_rings[i - 1][l0];
334 7554 : _gap_type[kgap] = EChannelType::CENTER;
335 7554 : kgap = kgap + 1;
336 7554 : _subch_type[k] = EChannelType::CENTER;
337 7554 : k = k + 1;
338 :
339 : } // for j
340 :
341 : // find the closest Pin at front ring
342 :
343 8097 : for (unsigned int j = 0; j < _pins_in_rings[i].size(); j++)
344 : {
345 7554 : if (j == _pins_in_rings[i].size() - 1)
346 : {
347 543 : _chan_to_pin_map[k].push_back(_pins_in_rings[i][j]);
348 543 : _chan_to_pin_map[k].push_back(_pins_in_rings[i][0]);
349 543 : avg_coor_x =
350 543 : 0.5 * (_pin_position[_pins_in_rings[i][j]](0) + _pin_position[_pins_in_rings[i][0]](0));
351 543 : avg_coor_y =
352 543 : 0.5 * (_pin_position[_pins_in_rings[i][j]](1) + _pin_position[_pins_in_rings[i][0]](1));
353 : }
354 : else
355 : {
356 7011 : _chan_to_pin_map[k].push_back(_pins_in_rings[i][j]);
357 7011 : _chan_to_pin_map[k].push_back(_pins_in_rings[i][j + 1]);
358 7011 : avg_coor_x = 0.5 * (_pin_position[_pins_in_rings[i][j]](0) +
359 7011 : _pin_position[_pins_in_rings[i][j + 1]](0));
360 7011 : avg_coor_y = 0.5 * (_pin_position[_pins_in_rings[i][j]](1) +
361 7011 : _pin_position[_pins_in_rings[i][j + 1]](1));
362 : }
363 :
364 : // if the outermost ring, set the edge subchannels first... then the corner subchannels
365 7554 : if (i == _n_rings - 1)
366 : {
367 : // add edges
368 3258 : _subch_type[k] = EChannelType::EDGE; // an edge subchannel is created
369 3258 : _gap_to_pin_map[kgap].first = _pins_in_rings[i][j];
370 3258 : _gap_to_pin_map[kgap].second = _pins_in_rings[i][j];
371 3258 : _gap_type[kgap] = EChannelType::EDGE;
372 3258 : _chan_to_gap_map[k].push_back(kgap);
373 3258 : kgap = kgap + 1;
374 3258 : k = k + 1;
375 :
376 3258 : if (j % i == 0)
377 : {
378 : // generate a corner subchannel, generate the additional gap and fix chan_to_gap_map
379 1104 : _gap_to_pin_map[kgap].first = _pins_in_rings[i][j];
380 1104 : _gap_to_pin_map[kgap].second = _pins_in_rings[i][j];
381 1104 : _gap_type[kgap] = EChannelType::CORNER;
382 :
383 : // corner subchannel
384 1104 : _chan_to_pin_map[k].push_back(_pins_in_rings[i][j]);
385 1104 : _chan_to_gap_map[k].push_back(kgap - 1);
386 1104 : _chan_to_gap_map[k].push_back(kgap);
387 1104 : _subch_type[k] = EChannelType::CORNER;
388 :
389 1104 : kgap = kgap + 1;
390 1104 : k = k + 1;
391 : }
392 : // if not the outer most ring
393 : }
394 : else
395 : {
396 : dist0 = 1.0e+5;
397 : unsigned int l0 = 0;
398 4296 : _chan_to_pin_map[k].push_back(_pins_in_rings[i + 1][0]);
399 105456 : for (unsigned int l = 0; l < _pins_in_rings[i + 1].size(); l++)
400 : {
401 101160 : dist = std::sqrt(pow(_pin_position[_pins_in_rings[i + 1][l]](0) - avg_coor_x, 2) +
402 101160 : pow(_pin_position[_pins_in_rings[i + 1][l]](1) - avg_coor_y, 2));
403 101160 : if (dist < dist0)
404 : {
405 31988 : _chan_to_pin_map[k][2] = _pins_in_rings[i + 1][l];
406 : dist0 = dist;
407 : l0 = l;
408 : } // if
409 : } // l
410 :
411 4296 : _gap_to_pin_map[kgap].first = _pins_in_rings[i][j];
412 4296 : _gap_to_pin_map[kgap].second = _pins_in_rings[i + 1][l0];
413 4296 : _gap_type[kgap] = EChannelType::CENTER;
414 4296 : kgap = kgap + 1;
415 4296 : _subch_type[k] = EChannelType::CENTER;
416 4296 : k = k + 1;
417 : } // if
418 : } // for j
419 : } // for i
420 :
421 : // Constructing pins to channels mao
422 7922 : for (unsigned int loc_rod = 0; loc_rod < _npins; loc_rod++)
423 : {
424 1313794 : for (unsigned int i = 0; i < _n_channels; i++)
425 : {
426 : bool rod_in_sc = false;
427 4931046 : for (unsigned int j : _chan_to_pin_map[i])
428 : {
429 3624990 : if (j == loc_rod)
430 : rod_in_sc = true;
431 : }
432 1306056 : if (rod_in_sc)
433 : {
434 43170 : _pin_to_chan_map[loc_rod].push_back(i);
435 : }
436 : }
437 : }
438 :
439 : // find the _gap_to_chan_map and _chan_to_gap_map using the gap_to_rod and subchannel_to_rod_maps
440 :
441 16396 : for (unsigned int i = 0; i < _n_channels; i++)
442 : {
443 16212 : if (_subch_type[i] == EChannelType::CENTER)
444 : {
445 3225318 : for (unsigned int j = 0; j < _n_gaps; j++)
446 : {
447 3213468 : if (_gap_type[j] == EChannelType::CENTER)
448 : {
449 2851632 : if (((_chan_to_pin_map[i][0] == _gap_to_pin_map[j].first) &&
450 2819340 : (_chan_to_pin_map[i][1] == _gap_to_pin_map[j].second)) ||
451 2808392 : ((_chan_to_pin_map[i][0] == _gap_to_pin_map[j].second) &&
452 31188 : (_chan_to_pin_map[i][1] == _gap_to_pin_map[j].first)))
453 : {
454 11850 : _chan_to_gap_map[i].push_back(j);
455 : }
456 :
457 2851632 : if (((_chan_to_pin_map[i][0] == _gap_to_pin_map[j].first) &&
458 2819340 : (_chan_to_pin_map[i][2] == _gap_to_pin_map[j].second)) ||
459 2807490 : ((_chan_to_pin_map[i][0] == _gap_to_pin_map[j].second) &&
460 31188 : (_chan_to_pin_map[i][2] == _gap_to_pin_map[j].first)))
461 : {
462 11850 : _chan_to_gap_map[i].push_back(j);
463 : }
464 :
465 2851632 : if (((_chan_to_pin_map[i][1] == _gap_to_pin_map[j].first) &&
466 2819340 : (_chan_to_pin_map[i][2] == _gap_to_pin_map[j].second)) ||
467 2816082 : ((_chan_to_pin_map[i][1] == _gap_to_pin_map[j].second) &&
468 31188 : (_chan_to_pin_map[i][2] == _gap_to_pin_map[j].first)))
469 : {
470 11850 : _chan_to_gap_map[i].push_back(j);
471 : }
472 : }
473 : } // for j
474 : }
475 4362 : else if (_subch_type[i] == EChannelType::EDGE)
476 : {
477 613998 : for (unsigned int j = 0; j < _n_gaps; j++)
478 : {
479 610740 : if (_gap_type[j] == EChannelType::CENTER)
480 : {
481 526608 : if (((_chan_to_pin_map[i][0] == _gap_to_pin_map[j].first) &&
482 520092 : (_chan_to_pin_map[i][1] == _gap_to_pin_map[j].second)) ||
483 517018 : ((_chan_to_pin_map[i][0] == _gap_to_pin_map[j].second) &&
484 5412 : (_chan_to_pin_map[i][1] == _gap_to_pin_map[j].first)))
485 : {
486 3258 : _chan_to_gap_map[i].push_back(j);
487 : }
488 : }
489 : }
490 :
491 : icorner = 0;
492 402771 : for (unsigned int k = 0; k < _n_channels; k++)
493 : {
494 400617 : if (_subch_type[k] == EChannelType::CORNER &&
495 16788 : _chan_to_pin_map[i][1] == _chan_to_pin_map[k][0])
496 : {
497 1104 : _chan_to_gap_map[i].push_back(_chan_to_gap_map[k][1]);
498 : icorner = 1;
499 : break;
500 : } // if
501 : } // for
502 :
503 402771 : for (unsigned int k = 0; k < _n_channels; k++)
504 : {
505 400617 : if (_subch_type[k] == EChannelType::CORNER &&
506 16788 : _chan_to_pin_map[i][0] == _chan_to_pin_map[k][0])
507 : {
508 1104 : _chan_to_gap_map[i].push_back(_chan_to_gap_map[k][1] + 1);
509 : icorner = 1;
510 : break;
511 : }
512 : }
513 :
514 2154 : if (icorner == 0)
515 : {
516 1050 : _chan_to_gap_map[i].push_back(_chan_to_gap_map[i][0] + 1);
517 : }
518 : }
519 : }
520 :
521 : // find gap_to_chan_map pair
522 :
523 23950 : for (unsigned int j = 0; j < _n_gaps; j++)
524 : {
525 3990570 : for (unsigned int i = 0; i < _n_channels; i++)
526 : {
527 3966804 : if (_subch_type[i] == EChannelType::CENTER || _subch_type[i] == EChannelType::EDGE)
528 : {
529 3824208 : if ((j == _chan_to_gap_map[i][0]) || (j == _chan_to_gap_map[i][1]) ||
530 3793992 : (j == _chan_to_gap_map[i][2]))
531 : {
532 45324 : if (_gap_to_chan_map[j].first == 0 && gap_fill[j].first == 0)
533 : {
534 23582 : _gap_to_chan_map[j].first = i;
535 23582 : gap_fill[j].first = 1;
536 : }
537 21742 : else if (_gap_to_chan_map[j].second == 0 && gap_fill[j].second == 0)
538 : {
539 21742 : _gap_to_chan_map[j].second = i;
540 21742 : gap_fill[j].second = 1;
541 : }
542 : else
543 : {
544 : }
545 : }
546 : }
547 142596 : else if (_subch_type[i] == EChannelType::CORNER)
548 : {
549 142596 : if ((j == _chan_to_gap_map[i][0]) || (j == _chan_to_gap_map[i][1]))
550 : {
551 2208 : if (_gap_to_chan_map[j].first == 0 && gap_fill[j].first == 0)
552 : {
553 184 : _gap_to_chan_map[j].first = i;
554 184 : gap_fill[j].first = 1;
555 : }
556 2024 : else if (_gap_to_chan_map[j].second == 0 && gap_fill[j].second == 0)
557 : {
558 2024 : _gap_to_chan_map[j].second = i;
559 2024 : gap_fill[j].second = 1;
560 : }
561 : else
562 : {
563 : }
564 : }
565 : }
566 : } // i
567 : } // j
568 :
569 16396 : for (unsigned int k = 0; k < _n_channels; k++)
570 : {
571 16212 : if (_subch_type[k] == EChannelType::EDGE)
572 : {
573 3258 : _gap_pairs_sf[k].first = _chan_to_gap_map[k][0];
574 3258 : _gap_pairs_sf[k].second = _chan_to_gap_map[k][2];
575 : auto k1 = _gap_pairs_sf[k].first;
576 : auto k2 = _gap_pairs_sf[k].second;
577 3258 : if (_gap_to_chan_map[k1].first == k)
578 : {
579 1104 : _chan_pairs_sf[k].first = _gap_to_chan_map[k1].second;
580 : }
581 : else
582 : {
583 2154 : _chan_pairs_sf[k].first = _gap_to_chan_map[k1].first;
584 : }
585 :
586 3258 : if (_gap_to_chan_map[k2].first == k)
587 : {
588 3074 : _chan_pairs_sf[k].second = _gap_to_chan_map[k2].second;
589 : }
590 : else
591 : {
592 184 : _chan_pairs_sf[k].second = _gap_to_chan_map[k2].first;
593 : }
594 : }
595 12954 : else if (_subch_type[k] == EChannelType::CORNER)
596 : {
597 1104 : _gap_pairs_sf[k].first = _chan_to_gap_map[k][1];
598 1104 : _gap_pairs_sf[k].second = _chan_to_gap_map[k][0];
599 :
600 : auto k1 = _gap_pairs_sf[k].first;
601 : auto k2 = _gap_pairs_sf[k].second;
602 :
603 1104 : if (_gap_to_chan_map[k1].first == k)
604 : {
605 184 : _chan_pairs_sf[k].first = _gap_to_chan_map[k1].second;
606 : }
607 : else
608 : {
609 920 : _chan_pairs_sf[k].first = _gap_to_chan_map[k1].first;
610 : }
611 :
612 1104 : if (_gap_to_chan_map[k2].first == k)
613 : {
614 0 : _chan_pairs_sf[k].second = _gap_to_chan_map[k2].second;
615 : }
616 : else
617 : {
618 1104 : _chan_pairs_sf[k].second = _gap_to_chan_map[k2].first;
619 : }
620 : }
621 : }
622 :
623 : // set the _gij_map
624 4705 : for (unsigned int iz = 0; iz < _n_cells + 1; iz++)
625 : {
626 593067 : for (unsigned int i_gap = 0; i_gap < _n_gaps; i_gap++)
627 : {
628 588546 : if (_gap_type[i_gap] == EChannelType::CENTER)
629 : {
630 481818 : _gij_map[iz].push_back(_pitch - _pin_diameter);
631 : }
632 106728 : else if (_gap_type[i_gap] == EChannelType::EDGE || _gap_type[i_gap] == EChannelType::CORNER)
633 : {
634 106728 : _gij_map[iz].push_back(_duct_to_pin_gap);
635 : }
636 : }
637 : }
638 :
639 16396 : for (unsigned int i = 0; i < _n_channels; i++)
640 : {
641 16212 : if (_subch_type[i] == EChannelType::CENTER || _subch_type[i] == EChannelType::EDGE)
642 : {
643 60432 : for (unsigned int k = 0; k < 3; k++)
644 : {
645 11517948 : for (unsigned int j = 0; j < _n_gaps; j++)
646 : {
647 11472624 : if (_chan_to_gap_map[i][k] == j && i == _gap_to_chan_map[j].first)
648 : {
649 23582 : if (i > _gap_to_chan_map[j].second)
650 : {
651 0 : _sign_id_crossflow_map[i][k] = negative_flow;
652 : }
653 : else
654 : {
655 23582 : _sign_id_crossflow_map[i][k] = positive_flow;
656 : }
657 : }
658 11449042 : else if (_chan_to_gap_map[i][k] == j && i == _gap_to_chan_map[j].second)
659 : {
660 21742 : if (i > _gap_to_chan_map[j].first)
661 : {
662 21742 : _sign_id_crossflow_map[i][k] = negative_flow;
663 : }
664 : else
665 : {
666 0 : _sign_id_crossflow_map[i][k] = positive_flow;
667 : }
668 : }
669 : } // j
670 : } // k
671 : }
672 1104 : else if (_subch_type[i] == EChannelType::CORNER)
673 : {
674 3312 : for (unsigned int k = 0; k < 2; k++)
675 : {
676 287400 : for (unsigned int j = 0; j < _n_gaps; j++)
677 : {
678 285192 : if (_chan_to_gap_map[i][k] == j && i == _gap_to_chan_map[j].first)
679 : {
680 184 : if (i > _gap_to_chan_map[j].second)
681 : {
682 0 : _sign_id_crossflow_map[i][k] = negative_flow;
683 : }
684 : else
685 : {
686 184 : _sign_id_crossflow_map[i][k] = positive_flow;
687 : }
688 : }
689 285008 : else if (_chan_to_gap_map[i][k] == j && i == _gap_to_chan_map[j].second)
690 : {
691 2024 : if (i > _gap_to_chan_map[j].first)
692 : {
693 2024 : _sign_id_crossflow_map[i][k] = negative_flow;
694 : }
695 : else
696 : {
697 0 : _sign_id_crossflow_map[i][k] = positive_flow;
698 : }
699 : }
700 : } // j
701 : } // k
702 : } // subch_type =2
703 : } // i
704 :
705 : // set the subchannel positions
706 16396 : for (unsigned int i = 0; i < _n_channels; i++)
707 : {
708 16212 : if (_subch_type[i] == EChannelType::CENTER)
709 : {
710 11850 : _subchannel_position[i][0] =
711 11850 : (_pin_position[_chan_to_pin_map[i][0]](0) + _pin_position[_chan_to_pin_map[i][1]](0) +
712 11850 : _pin_position[_chan_to_pin_map[i][2]](0)) /
713 : 3.0;
714 11850 : _subchannel_position[i][1] =
715 11850 : (_pin_position[_chan_to_pin_map[i][0]](1) + _pin_position[_chan_to_pin_map[i][1]](1) +
716 11850 : _pin_position[_chan_to_pin_map[i][2]](1)) /
717 : 3.0;
718 : }
719 4362 : else if (_subch_type[i] == EChannelType::EDGE)
720 : {
721 416934 : for (unsigned int j = 0; j < _n_channels; j++)
722 : {
723 413676 : if (_subch_type[j] == EChannelType::CENTER &&
724 323028 : ((_chan_to_pin_map[i][0] == _chan_to_pin_map[j][0] &&
725 3258 : _chan_to_pin_map[i][1] == _chan_to_pin_map[j][1]) ||
726 319770 : (_chan_to_pin_map[i][0] == _chan_to_pin_map[j][1] &&
727 3258 : _chan_to_pin_map[i][1] == _chan_to_pin_map[j][0])))
728 : {
729 3258 : x0 = _pin_position[_chan_to_pin_map[j][2]](0);
730 3258 : y0 = _pin_position[_chan_to_pin_map[j][2]](1);
731 : }
732 410418 : else if (_subch_type[j] == EChannelType::CENTER &&
733 319770 : ((_chan_to_pin_map[i][0] == _chan_to_pin_map[j][0] &&
734 0 : _chan_to_pin_map[i][1] == _chan_to_pin_map[j][2]) ||
735 319770 : (_chan_to_pin_map[i][0] == _chan_to_pin_map[j][2] &&
736 2154 : _chan_to_pin_map[i][1] == _chan_to_pin_map[j][0])))
737 : {
738 0 : x0 = _pin_position[_chan_to_pin_map[j][1]](0);
739 0 : y0 = _pin_position[_chan_to_pin_map[j][1]](1);
740 : }
741 410418 : else if (_subch_type[j] == EChannelType::CENTER &&
742 319770 : ((_chan_to_pin_map[i][0] == _chan_to_pin_map[j][1] &&
743 3258 : _chan_to_pin_map[i][1] == _chan_to_pin_map[j][2]) ||
744 319770 : (_chan_to_pin_map[i][0] == _chan_to_pin_map[j][2] &&
745 2154 : _chan_to_pin_map[i][1] == _chan_to_pin_map[j][1])))
746 : {
747 0 : x0 = _pin_position[_chan_to_pin_map[j][0]](0);
748 0 : y0 = _pin_position[_chan_to_pin_map[j][0]](1);
749 : }
750 413676 : x1 = 0.5 *
751 413676 : (_pin_position[_chan_to_pin_map[i][0]](0) + _pin_position[_chan_to_pin_map[i][1]](0));
752 413676 : y1 = 0.5 *
753 413676 : (_pin_position[_chan_to_pin_map[i][0]](1) + _pin_position[_chan_to_pin_map[i][1]](1));
754 413676 : a1 = _pin_diameter / 2.0 + _duct_to_pin_gap / 2.0;
755 413676 : a2 = std::sqrt((x1 - x0) * (x1 - x0) + (y1 - y0) * (y1 - y0)) + a1;
756 413676 : _subchannel_position[i][0] = (a2 * x1 - a1 * x0) / (a2 - a1);
757 413676 : _subchannel_position[i][1] = (a2 * y1 - a1 * y0) / (a2 - a1);
758 : } // j
759 : }
760 1104 : else if (_subch_type[i] == EChannelType::CORNER)
761 : {
762 1104 : x0 = _pin_position[0](0);
763 1104 : y0 = _pin_position[0](1);
764 1104 : x1 = _pin_position[_chan_to_pin_map[i][0]](0);
765 1104 : y1 = _pin_position[_chan_to_pin_map[i][0]](1);
766 1104 : a1 = _pin_diameter / 2.0 + _duct_to_pin_gap / 2.0;
767 1104 : a2 = std::sqrt((x1 - x0) * (x1 - x0) + (y1 - y0) * (y1 - y0)) + a1;
768 1104 : _subchannel_position[i][0] = (a2 * x1 - a1 * x0) / (a2 - a1);
769 1104 : _subchannel_position[i][1] = (a2 * y1 - a1 * y0) / (a2 - a1);
770 : }
771 : }
772 :
773 : // Reduce reserved memory in the channel-to-gap map.
774 16396 : for (auto & gap : _chan_to_gap_map)
775 : {
776 : gap.shrink_to_fit();
777 : }
778 184 : }
779 :
780 : void
781 184 : SCMTriAssemblyMeshGenerator::buildPinMesh(MeshBase & mesh_base)
782 : {
783 184 : if (_npins == 0)
784 : return;
785 :
786 184 : mesh_base.reserve_elem(_n_cells * _npins);
787 184 : mesh_base.reserve_nodes((_n_cells + 1) * _npins);
788 :
789 184 : _pin_nodes.clear();
790 184 : _pin_nodes.resize(_npins);
791 :
792 : // Defining the extent of the subchannel mesh to append the pin mesh to the current subchannel
793 : // mesh.
794 184 : const unsigned int node_sub = mesh_base.n_nodes();
795 184 : const unsigned int elem_sub = mesh_base.n_elem();
796 :
797 : // Add the points in the shape of a rectilinear grid. The grid is regular on the xy-plane at the
798 : // triangular lattice pin positions. The grid along z is also regular. Store pointers in the
799 : // _pin_nodes array so we can keep track of which points are in which pins.
800 : unsigned int node_id = node_sub;
801 7922 : for (unsigned int i = 0; i < _npins; i++)
802 : {
803 7738 : _pin_nodes[i].reserve(_n_cells + 1);
804 199399 : for (unsigned int iz = 0; iz < _n_cells + 1; iz++)
805 191661 : _pin_nodes[i].push_back(mesh_base.add_point(
806 383322 : Point(_pin_position[i](0), _pin_position[i](1), _z_grid[iz]), node_id++));
807 : }
808 :
809 : // Add the elements which in this case are 2-node edges that link each pin's nodes vertically.
810 : unsigned int elem_id = elem_sub;
811 7922 : for (unsigned int i = 0; i < _npins; i++)
812 191661 : for (unsigned int iz = 0; iz < _n_cells; iz++)
813 : {
814 183923 : Elem * elem = mesh_base.add_elem(std::make_unique<Edge2>());
815 183923 : elem->subdomain_id() = _pin_block_id;
816 183923 : elem->set_id(elem_id++);
817 183923 : const int indx1 = (_n_cells + 1) * i + iz + node_sub;
818 183923 : const int indx2 = (_n_cells + 1) * i + (iz + 1) + node_sub;
819 183923 : elem->set_node(0, mesh_base.node_ptr(indx1));
820 183923 : elem->set_node(1, mesh_base.node_ptr(indx2));
821 : }
822 :
823 184 : mesh_base.subdomain_name(_pin_block_id) = "fuel_pins";
824 : }
825 :
826 : std::unique_ptr<MeshBase>
827 184 : SCMTriAssemblyMeshGenerator::generate()
828 : {
829 184 : auto mesh_base = buildMeshBaseObject();
830 :
831 : BoundaryInfo & boundary_info = mesh_base->get_boundary_info();
832 184 : mesh_base->set_spatial_dimension(3);
833 184 : mesh_base->reserve_elem(_n_cells * (_n_channels + _npins));
834 184 : mesh_base->reserve_nodes((_n_cells + 1) * (_n_channels + _npins));
835 184 : _nodes.resize(_n_channels);
836 : // Add the points for the give x,y subchannel positions. The grid is hexagonal.
837 : // The grid along
838 : // z is irregular to account for Pin spacers. Store pointers in the _nodes
839 : // array so we can keep track of which points are in which channels.
840 : unsigned int node_id = 0;
841 16396 : for (unsigned int i = 0; i < _n_channels; i++)
842 : {
843 16212 : _nodes[i].reserve(_n_cells + 1);
844 417618 : for (unsigned int iz = 0; iz < _n_cells + 1; iz++)
845 : {
846 401406 : _nodes[i].push_back(mesh_base->add_point(
847 802812 : Point(_subchannel_position[i][0], _subchannel_position[i][1], _z_grid[iz]), node_id++));
848 : }
849 : }
850 :
851 : // Add the elements which in this case are 2-node edges that link each
852 : // subchannel's nodes vertically.
853 : unsigned int elem_id = 0;
854 16396 : for (unsigned int i = 0; i < _n_channels; i++)
855 : {
856 401406 : for (unsigned int iz = 0; iz < _n_cells; iz++)
857 : {
858 385194 : Elem * elem = mesh_base->add_elem(std::make_unique<Edge2>());
859 385194 : elem->subdomain_id() = _subchannel_block_id;
860 385194 : elem->set_id(elem_id++);
861 385194 : const int indx1 = (_n_cells + 1) * i + iz;
862 385194 : const int indx2 = (_n_cells + 1) * i + (iz + 1);
863 385194 : elem->set_node(0, mesh_base->node_ptr(indx1));
864 385194 : elem->set_node(1, mesh_base->node_ptr(indx2));
865 :
866 385194 : if (iz == 0)
867 16212 : boundary_info.add_side(elem, 0, 0);
868 385194 : if (iz == _n_cells - 1)
869 16212 : boundary_info.add_side(elem, 1, 1);
870 : }
871 : }
872 184 : boundary_info.sideset_name(0) = "inlet";
873 184 : boundary_info.sideset_name(1) = "outlet";
874 184 : boundary_info.nodeset_name(0) = "inlet";
875 184 : boundary_info.nodeset_name(1) = "outlet";
876 :
877 : // Naming the block
878 184 : mesh_base->subdomain_name(_subchannel_block_id) = "subchannel";
879 184 : buildPinMesh(*mesh_base);
880 :
881 184 : mesh_base->prepare_for_use();
882 :
883 : // move the meta data into TriSubChannelMesh
884 184 : auto & sch_mesh = static_cast<TriSubChannelMesh &>(*_mesh);
885 184 : sch_mesh._unheated_length_entry = _unheated_length_entry;
886 184 : sch_mesh._heated_length = _heated_length;
887 184 : sch_mesh._unheated_length_exit = _unheated_length_exit;
888 184 : sch_mesh._z_grid = _z_grid;
889 184 : sch_mesh._k_grid = _k_grid;
890 184 : sch_mesh._spacer_z = _spacer_z;
891 184 : sch_mesh._spacer_k = _spacer_k;
892 184 : sch_mesh._z_blockage = _z_blockage;
893 184 : sch_mesh._index_blockage = _index_blockage;
894 184 : sch_mesh._reduction_blockage = _reduction_blockage;
895 184 : sch_mesh._kij = _kij;
896 184 : sch_mesh._pitch = _pitch;
897 184 : sch_mesh._pin_diameter = _pin_diameter;
898 184 : sch_mesh._n_cells = _n_cells;
899 184 : sch_mesh._n_rings = _n_rings;
900 184 : sch_mesh._n_channels = _n_channels;
901 184 : sch_mesh._flat_to_flat = _flat_to_flat;
902 184 : sch_mesh._dwire = _dwire;
903 184 : sch_mesh._hwire = _hwire;
904 184 : sch_mesh._duct_to_pin_gap = _duct_to_pin_gap;
905 184 : sch_mesh._nodes = _nodes;
906 184 : sch_mesh._gap_to_chan_map = _gap_to_chan_map;
907 184 : sch_mesh._gap_to_pin_map = _gap_to_pin_map;
908 184 : sch_mesh._chan_to_gap_map = _chan_to_gap_map;
909 184 : sch_mesh._sign_id_crossflow_map = _sign_id_crossflow_map;
910 184 : sch_mesh._gij_map = _gij_map;
911 184 : sch_mesh._subchannel_position = _subchannel_position;
912 184 : sch_mesh._pin_position = _pin_position;
913 184 : sch_mesh._pins_in_rings = _pins_in_rings;
914 184 : sch_mesh._chan_to_pin_map = _chan_to_pin_map;
915 184 : sch_mesh._npins = _npins;
916 184 : sch_mesh._n_gaps = _n_gaps;
917 184 : sch_mesh._subch_type = _subch_type;
918 184 : sch_mesh._gap_type = _gap_type;
919 184 : sch_mesh._gap_pairs_sf = _gap_pairs_sf;
920 184 : sch_mesh._chan_pairs_sf = _chan_pairs_sf;
921 184 : sch_mesh._pin_to_chan_map = _pin_to_chan_map;
922 184 : sch_mesh._pin_nodes = _pin_nodes;
923 184 : sch_mesh._pin_mesh_exist = (_npins > 0);
924 184 : sch_mesh.computeAssemblyHydraulicParameters();
925 :
926 184 : return mesh_base;
927 0 : }
|
