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