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 332 : SCMTriSubChannelMeshGenerator::validParams()
24 : {
25 332 : InputParameters params = MeshGenerator::validParams();
26 332 : params.addClassDescription(
27 : "Creates a mesh of 1D subchannels in a triangular lattice arrangement");
28 664 : params.addRequiredParam<unsigned int>("n_cells", "The number of cells in the axial direction");
29 664 : params.addRequiredParam<Real>("pitch", "Pitch [m]");
30 664 : params.addRequiredParam<Real>("pin_diameter", "Rod diameter [m]");
31 664 : params.addParam<Real>("unheated_length_entry", 0.0, "Unheated length at entry [m]");
32 664 : params.addRequiredParam<Real>("heated_length", "Heated length [m]");
33 664 : params.addParam<Real>("unheated_length_exit", 0.0, "Unheated length at exit [m]");
34 664 : params.addRequiredParam<unsigned int>("nrings", "Number of fuel Pin rings per assembly [-]");
35 664 : params.addRequiredParam<Real>("flat_to_flat",
36 : "Flat to flat distance for the hexagonal assembly [m]");
37 664 : params.addRequiredParam<Real>("dwire", "Wire diameter [m]");
38 664 : params.addRequiredParam<Real>("hwire", "Wire lead length [m]");
39 664 : params.addParam<std::vector<Real>>(
40 : "spacer_z", {}, "Axial location of spacers/vanes/mixing_vanes [m]");
41 664 : params.addParam<std::vector<Real>>(
42 : "spacer_k", {}, "K-loss coefficient of spacers/vanes/mixing_vanes [-]");
43 664 : params.addParam<Real>("Kij", 0.5, "Lateral form loss coefficient [-]");
44 664 : params.addParam<std::vector<Real>>("z_blockage",
45 664 : std::vector<Real>({0.0, 0.0}),
46 : "axial location of blockage (inlet, outlet) [m]");
47 664 : params.addParam<std::vector<unsigned int>>("index_blockage",
48 664 : std::vector<unsigned int>({0}),
49 : "index of subchannels affected by blockage");
50 664 : params.addParam<std::vector<Real>>(
51 : "reduction_blockage",
52 664 : std::vector<Real>({1.0}),
53 : "Area reduction of subchannels affected by blockage (number to muliply the area)");
54 664 : params.addParam<std::vector<Real>>("k_blockage",
55 664 : std::vector<Real>({0.0}),
56 : "Form loss coefficient of subchannels affected by blockage");
57 664 : params.addParam<unsigned int>("block_id", 0, "Domain Index");
58 332 : return params;
59 0 : }
60 :
61 166 : SCMTriSubChannelMeshGenerator::SCMTriSubChannelMeshGenerator(const InputParameters & params)
62 : : MeshGenerator(params),
63 166 : _unheated_length_entry(getParam<Real>("unheated_length_entry")),
64 332 : _heated_length(getParam<Real>("heated_length")),
65 332 : _unheated_length_exit(getParam<Real>("unheated_length_exit")),
66 332 : _block_id(getParam<unsigned int>("block_id")),
67 332 : _spacer_z(getParam<std::vector<Real>>("spacer_z")),
68 332 : _spacer_k(getParam<std::vector<Real>>("spacer_k")),
69 332 : _z_blockage(getParam<std::vector<Real>>("z_blockage")),
70 332 : _index_blockage(getParam<std::vector<unsigned int>>("index_blockage")),
71 332 : _reduction_blockage(getParam<std::vector<Real>>("reduction_blockage")),
72 332 : _k_blockage(getParam<std::vector<Real>>("k_blockage")),
73 332 : _pitch(getParam<Real>("pitch")),
74 332 : _kij(getParam<Real>("Kij")),
75 332 : _pin_diameter(getParam<Real>("pin_diameter")),
76 332 : _n_cells(getParam<unsigned int>("n_cells")),
77 332 : _n_rings(getParam<unsigned int>("nrings")),
78 332 : _flat_to_flat(getParam<Real>("flat_to_flat")),
79 332 : _dwire(getParam<Real>("dwire")),
80 332 : _hwire(getParam<Real>("hwire")),
81 166 : _duct_to_pin_gap(0.5 *
82 332 : (_flat_to_flat - (_n_rings - 1) * _pitch * std::sqrt(3.0) - _pin_diameter))
83 : {
84 166 : 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 166 : if (_spacer_z.size() &&
88 124 : _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 166 : if (_z_blockage.size() != 2)
92 0 : mooseError(name(), ": Size of vector z_blockage must be 2");
93 :
94 166 : 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 166 : if ((_index_blockage.size() != _reduction_blockage.size()) ||
98 332 : (_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 166 : SubChannelMesh::generateZGrid(
105 166 : _unheated_length_entry, _heated_length, _unheated_length_exit, _n_cells, _z_grid);
106 :
107 : // Defining the total length from 3 axial sections
108 166 : 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 338 : for (const auto & elem : _spacer_z)
113 172 : spacer_cell.emplace_back(std::round(elem * _n_cells / L));
114 :
115 : // Defining the array for axial resistances
116 : std::vector<Real> kgrid;
117 166 : kgrid.resize(_n_cells + 1, 0.0);
118 :
119 : // Summing the spacer resistance to the grid resistance array
120 338 : for (unsigned int index = 0; index < spacer_cell.size(); index++)
121 172 : 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 166 : 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 166 : 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 166 : TriSubChannelMesh::rodPositions(_pin_position, _n_rings, _pitch, Point(0, 0));
157 166 : _npins = _pin_position.size();
158 : // assign the pins to the corresponding rings
159 : unsigned int k = 0; // initialize the fuel Pin counter index
160 166 : _pins_in_rings.resize(_n_rings);
161 166 : _pins_in_rings[0].push_back(k++);
162 583 : for (unsigned int i = 1; i < _n_rings; i++)
163 5415 : for (unsigned int j = 0; j < i * 6; j++)
164 4998 : _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 583 : for (unsigned int j = 0; j < _n_rings - 1; j++)
170 417 : chancount += j * 6;
171 : // Adding external channels to the total count
172 166 : _n_channels = chancount + _npins - 1 + (_n_rings - 1) * 6 + 6;
173 :
174 166 : 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 166 : 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 166 : _k_grid.resize(_n_channels, std::vector<Real>(_n_cells + 1));
187 11158 : for (unsigned int i = 0; i < _n_channels; i++)
188 10992 : _k_grid[i] = kgrid;
189 :
190 : // Add blockage resistance to the 2D grid resistane array
191 166 : Real dz = L / _n_cells;
192 4144 : for (unsigned int i = 0; i < _n_cells + 1; i++)
193 : {
194 3978 : if ((dz * i >= _z_blockage.front() && dz * i <= _z_blockage.back()))
195 : {
196 : unsigned int index(0);
197 518 : for (const auto & i_ch : _index_blockage)
198 : {
199 352 : _k_grid[i_ch][i] += _k_blockage[index];
200 352 : index++;
201 : }
202 : }
203 : }
204 :
205 166 : _chan_to_pin_map.resize(_n_channels);
206 166 : _pin_to_chan_map.resize(_npins);
207 166 : _subch_type.resize(_n_channels);
208 166 : _n_gaps = _n_channels + _npins - 1; /// initial assignment
209 166 : _gap_to_chan_map.resize(_n_gaps);
210 166 : _gap_to_pin_map.resize(_n_gaps);
211 166 : gap_fill.resize(_n_gaps);
212 166 : _chan_to_gap_map.resize(_n_channels);
213 166 : _gap_pairs_sf.resize(_n_channels);
214 166 : _chan_pairs_sf.resize(_n_channels);
215 166 : _gij_map.resize(_n_cells + 1);
216 166 : _sign_id_crossflow_map.resize(_n_channels);
217 166 : _gap_type.resize(_n_gaps);
218 166 : _subchannel_position.resize(_n_channels);
219 :
220 11158 : for (unsigned int i = 0; i < _n_channels; i++)
221 : {
222 10992 : _chan_to_pin_map[i].reserve(3);
223 10992 : _chan_to_gap_map[i].reserve(3);
224 10992 : _sign_id_crossflow_map[i].reserve(3);
225 10992 : _subchannel_position[i].reserve(3);
226 43968 : for (unsigned int j = 0; j < 3; j++)
227 : {
228 32976 : _sign_id_crossflow_map.at(i).push_back(positive_flow);
229 32976 : _subchannel_position.at(i).push_back(0.0);
230 : }
231 : }
232 :
233 4144 : for (unsigned int iz = 0; iz < _n_cells + 1; iz++)
234 : {
235 3978 : _gij_map[iz].reserve(_n_gaps);
236 : }
237 :
238 5330 : for (unsigned int i = 0; i < _npins; i++)
239 5164 : _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 583 : for (unsigned int i = 1; i < _n_rings; i++)
247 : {
248 : // find the closest Pin at back ring
249 5415 : for (unsigned int j = 0; j < _pins_in_rings[i].size(); j++)
250 : {
251 4998 : if (j == _pins_in_rings[i].size() - 1)
252 : {
253 417 : _chan_to_pin_map[k].push_back(_pins_in_rings[i][j]);
254 417 : _chan_to_pin_map[k].push_back(_pins_in_rings[i][0]);
255 417 : avg_coor_x =
256 417 : 0.5 * (_pin_position[_pins_in_rings[i][j]](0) + _pin_position[_pins_in_rings[i][0]](0));
257 417 : avg_coor_y =
258 417 : 0.5 * (_pin_position[_pins_in_rings[i][j]](1) + _pin_position[_pins_in_rings[i][0]](1));
259 417 : _gap_to_pin_map[kgap].first = _pins_in_rings[i][0];
260 417 : _gap_to_pin_map[kgap].second = _pins_in_rings[i][j];
261 417 : _gap_type[kgap] = EChannelType::CENTER;
262 417 : kgap = kgap + 1;
263 : }
264 : else
265 : {
266 4581 : _chan_to_pin_map[k].push_back(_pins_in_rings[i][j]);
267 4581 : _chan_to_pin_map[k].push_back(_pins_in_rings[i][j + 1]);
268 4581 : avg_coor_x = 0.5 * (_pin_position[_pins_in_rings[i][j]](0) +
269 4581 : _pin_position[_pins_in_rings[i][j + 1]](0));
270 4581 : avg_coor_y = 0.5 * (_pin_position[_pins_in_rings[i][j]](1) +
271 4581 : _pin_position[_pins_in_rings[i][j + 1]](1));
272 4581 : _gap_to_pin_map[kgap].first = _pins_in_rings[i][j];
273 4581 : _gap_to_pin_map[kgap].second = _pins_in_rings[i][j + 1];
274 4581 : _gap_type[kgap] = EChannelType::CENTER;
275 4581 : kgap = kgap + 1;
276 : }
277 :
278 : dist0 = 1.0e+5;
279 :
280 4998 : _chan_to_pin_map[k].push_back(_pins_in_rings[i - 1][0]);
281 : unsigned int l0 = 0;
282 :
283 57906 : for (unsigned int l = 0; l < _pins_in_rings[i - 1].size(); l++)
284 : {
285 52908 : dist = std::sqrt(pow(_pin_position[_pins_in_rings[i - 1][l]](0) - avg_coor_x, 2) +
286 52908 : pow(_pin_position[_pins_in_rings[i - 1][l]](1) - avg_coor_y, 2));
287 :
288 52908 : if (dist < dist0)
289 : {
290 18986 : _chan_to_pin_map[k][2] = _pins_in_rings[i - 1][l];
291 : l0 = l;
292 : dist0 = dist;
293 : } // if
294 : } // l
295 :
296 4998 : _gap_to_pin_map[kgap].first = _pins_in_rings[i][j];
297 4998 : _gap_to_pin_map[kgap].second = _pins_in_rings[i - 1][l0];
298 4998 : _gap_type[kgap] = EChannelType::CENTER;
299 4998 : kgap = kgap + 1;
300 4998 : _subch_type[k] = EChannelType::CENTER;
301 4998 : k = k + 1;
302 :
303 : } // for j
304 :
305 : // find the closest Pin at front ring
306 :
307 5415 : for (unsigned int j = 0; j < _pins_in_rings[i].size(); j++)
308 : {
309 4998 : if (j == _pins_in_rings[i].size() - 1)
310 : {
311 417 : _chan_to_pin_map[k].push_back(_pins_in_rings[i][j]);
312 417 : _chan_to_pin_map[k].push_back(_pins_in_rings[i][0]);
313 417 : avg_coor_x =
314 417 : 0.5 * (_pin_position[_pins_in_rings[i][j]](0) + _pin_position[_pins_in_rings[i][0]](0));
315 417 : avg_coor_y =
316 417 : 0.5 * (_pin_position[_pins_in_rings[i][j]](1) + _pin_position[_pins_in_rings[i][0]](1));
317 : }
318 : else
319 : {
320 4581 : _chan_to_pin_map[k].push_back(_pins_in_rings[i][j]);
321 4581 : _chan_to_pin_map[k].push_back(_pins_in_rings[i][j + 1]);
322 4581 : avg_coor_x = 0.5 * (_pin_position[_pins_in_rings[i][j]](0) +
323 4581 : _pin_position[_pins_in_rings[i][j + 1]](0));
324 4581 : avg_coor_y = 0.5 * (_pin_position[_pins_in_rings[i][j]](1) +
325 4581 : _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 4998 : if (i == _n_rings - 1)
330 : {
331 : // add edges
332 2502 : _subch_type[k] = EChannelType::EDGE; // an edge subchannel is created
333 2502 : _gap_to_pin_map[kgap].first = _pins_in_rings[i][j];
334 2502 : _gap_to_pin_map[kgap].second = _pins_in_rings[i][j];
335 2502 : _gap_type[kgap] = EChannelType::EDGE;
336 2502 : _chan_to_gap_map[k].push_back(kgap);
337 2502 : kgap = kgap + 1;
338 2502 : k = k + 1;
339 :
340 2502 : if (j % i == 0)
341 : {
342 : // generate a corner subchannel, generate the additional gap and fix chan_to_gap_map
343 996 : _gap_to_pin_map[kgap].first = _pins_in_rings[i][j];
344 996 : _gap_to_pin_map[kgap].second = _pins_in_rings[i][j];
345 996 : _gap_type[kgap] = EChannelType::CORNER;
346 :
347 : // corner subchannel
348 996 : _chan_to_pin_map[k].push_back(_pins_in_rings[i][j]);
349 996 : _chan_to_gap_map[k].push_back(kgap - 1);
350 996 : _chan_to_gap_map[k].push_back(kgap);
351 996 : _subch_type[k] = EChannelType::CORNER;
352 :
353 996 : kgap = kgap + 1;
354 996 : 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 2496 : _chan_to_pin_map[k].push_back(_pins_in_rings[i + 1][0]);
363 54408 : for (unsigned int l = 0; l < _pins_in_rings[i + 1].size(); l++)
364 : {
365 51912 : dist = std::sqrt(pow(_pin_position[_pins_in_rings[i + 1][l]](0) - avg_coor_x, 2) +
366 51912 : pow(_pin_position[_pins_in_rings[i + 1][l]](1) - avg_coor_y, 2));
367 51912 : if (dist < dist0)
368 : {
369 16801 : _chan_to_pin_map[k][2] = _pins_in_rings[i + 1][l];
370 : dist0 = dist;
371 : l0 = l;
372 : } // if
373 : } // l
374 :
375 2496 : _gap_to_pin_map[kgap].first = _pins_in_rings[i][j];
376 2496 : _gap_to_pin_map[kgap].second = _pins_in_rings[i + 1][l0];
377 2496 : _gap_type[kgap] = EChannelType::CENTER;
378 2496 : kgap = kgap + 1;
379 2496 : _subch_type[k] = EChannelType::CENTER;
380 2496 : k = k + 1;
381 : } // if
382 : } // for j
383 : } // for i
384 :
385 : // Constructing pins to channels mao
386 5330 : for (unsigned int loc_rod = 0; loc_rod < _npins; loc_rod++)
387 : {
388 648712 : for (unsigned int i = 0; i < _n_channels; i++)
389 : {
390 : bool rod_in_sc = false;
391 2401866 : for (unsigned int j : _chan_to_pin_map[i])
392 : {
393 1758318 : if (j == loc_rod)
394 : rod_in_sc = true;
395 : }
396 643548 : if (rod_in_sc)
397 : {
398 28482 : _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 11158 : for (unsigned int i = 0; i < _n_channels; i++)
406 : {
407 10992 : if (_subch_type[i] == EChannelType::CENTER)
408 : {
409 1536594 : for (unsigned int j = 0; j < _n_gaps; j++)
410 : {
411 1529100 : if (_gap_type[j] == EChannelType::CENTER)
412 : {
413 1333368 : if (((_chan_to_pin_map[i][0] == _gap_to_pin_map[j].first) &&
414 1313388 : (_chan_to_pin_map[i][1] == _gap_to_pin_map[j].second)) ||
415 1306562 : ((_chan_to_pin_map[i][0] == _gap_to_pin_map[j].second) &&
416 18984 : (_chan_to_pin_map[i][1] == _gap_to_pin_map[j].first)))
417 : {
418 7494 : _chan_to_gap_map[i].push_back(j);
419 : }
420 :
421 1333368 : if (((_chan_to_pin_map[i][0] == _gap_to_pin_map[j].first) &&
422 1313388 : (_chan_to_pin_map[i][2] == _gap_to_pin_map[j].second)) ||
423 1305894 : ((_chan_to_pin_map[i][0] == _gap_to_pin_map[j].second) &&
424 18984 : (_chan_to_pin_map[i][2] == _gap_to_pin_map[j].first)))
425 : {
426 7494 : _chan_to_gap_map[i].push_back(j);
427 : }
428 :
429 1333368 : if (((_chan_to_pin_map[i][1] == _gap_to_pin_map[j].first) &&
430 1313388 : (_chan_to_pin_map[i][2] == _gap_to_pin_map[j].second)) ||
431 1310886 : ((_chan_to_pin_map[i][1] == _gap_to_pin_map[j].second) &&
432 18984 : (_chan_to_pin_map[i][2] == _gap_to_pin_map[j].first)))
433 : {
434 7494 : _chan_to_gap_map[i].push_back(j);
435 : }
436 : }
437 : } // for j
438 : }
439 3498 : else if (_subch_type[i] == EChannelType::EDGE)
440 : {
441 341082 : for (unsigned int j = 0; j < _n_gaps; j++)
442 : {
443 338580 : if (_gap_type[j] == EChannelType::CENTER)
444 : {
445 283608 : if (((_chan_to_pin_map[i][0] == _gap_to_pin_map[j].first) &&
446 278604 : (_chan_to_pin_map[i][1] == _gap_to_pin_map[j].second)) ||
447 276268 : ((_chan_to_pin_map[i][0] == _gap_to_pin_map[j].second) &&
448 4008 : (_chan_to_pin_map[i][1] == _gap_to_pin_map[j].first)))
449 : {
450 2502 : _chan_to_gap_map[i].push_back(j);
451 : }
452 : }
453 : }
454 :
455 : icorner = 0;
456 221979 : for (unsigned int k = 0; k < _n_channels; k++)
457 : {
458 220473 : if (_subch_type[k] == EChannelType::CORNER &&
459 12522 : _chan_to_pin_map[i][1] == _chan_to_pin_map[k][0])
460 : {
461 996 : _chan_to_gap_map[i].push_back(_chan_to_gap_map[k][1]);
462 : icorner = 1;
463 : break;
464 : } // if
465 : } // for
466 :
467 221979 : for (unsigned int k = 0; k < _n_channels; k++)
468 : {
469 220473 : if (_subch_type[k] == EChannelType::CORNER &&
470 12522 : _chan_to_pin_map[i][0] == _chan_to_pin_map[k][0])
471 : {
472 996 : _chan_to_gap_map[i].push_back(_chan_to_gap_map[k][1] + 1);
473 : icorner = 1;
474 : break;
475 : }
476 : }
477 :
478 1506 : if (icorner == 0)
479 : {
480 510 : _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 16156 : for (unsigned int j = 0; j < _n_gaps; j++)
488 : {
489 1979610 : for (unsigned int i = 0; i < _n_channels; i++)
490 : {
491 1963620 : if (_subch_type[i] == EChannelType::CENTER || _subch_type[i] == EChannelType::EDGE)
492 : {
493 1867680 : if ((j == _chan_to_gap_map[i][0]) || (j == _chan_to_gap_map[i][1]) ||
494 1847688 : (j == _chan_to_gap_map[i][2]))
495 : {
496 29988 : if (_gap_to_chan_map[j].first == 0 && gap_fill[j].first == 0)
497 : {
498 15824 : _gap_to_chan_map[j].first = i;
499 15824 : gap_fill[j].first = 1;
500 : }
501 14164 : else if (_gap_to_chan_map[j].second == 0 && gap_fill[j].second == 0)
502 : {
503 14164 : _gap_to_chan_map[j].second = i;
504 14164 : gap_fill[j].second = 1;
505 : }
506 : else
507 : {
508 : }
509 : }
510 : }
511 95940 : else if (_subch_type[i] == EChannelType::CORNER)
512 : {
513 95940 : if ((j == _chan_to_gap_map[i][0]) || (j == _chan_to_gap_map[i][1]))
514 : {
515 1992 : if (_gap_to_chan_map[j].first == 0 && gap_fill[j].first == 0)
516 : {
517 166 : _gap_to_chan_map[j].first = i;
518 166 : gap_fill[j].first = 1;
519 : }
520 1826 : else if (_gap_to_chan_map[j].second == 0 && gap_fill[j].second == 0)
521 : {
522 1826 : _gap_to_chan_map[j].second = i;
523 1826 : gap_fill[j].second = 1;
524 : }
525 : else
526 : {
527 : }
528 : }
529 : }
530 : } // i
531 : } // j
532 :
533 11158 : for (unsigned int k = 0; k < _n_channels; k++)
534 : {
535 10992 : if (_subch_type[k] == EChannelType::EDGE)
536 : {
537 2502 : _gap_pairs_sf[k].first = _chan_to_gap_map[k][0];
538 2502 : _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 2502 : if (_gap_to_chan_map[k1].first == k)
542 : {
543 996 : _chan_pairs_sf[k].first = _gap_to_chan_map[k1].second;
544 : }
545 : else
546 : {
547 1506 : _chan_pairs_sf[k].first = _gap_to_chan_map[k1].first;
548 : }
549 :
550 2502 : if (_gap_to_chan_map[k2].first == k)
551 : {
552 2336 : _chan_pairs_sf[k].second = _gap_to_chan_map[k2].second;
553 : }
554 : else
555 : {
556 166 : _chan_pairs_sf[k].second = _gap_to_chan_map[k2].first;
557 : }
558 : }
559 8490 : else if (_subch_type[k] == EChannelType::CORNER)
560 : {
561 996 : _gap_pairs_sf[k].first = _chan_to_gap_map[k][1];
562 996 : _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 996 : if (_gap_to_chan_map[k1].first == k)
568 : {
569 166 : _chan_pairs_sf[k].first = _gap_to_chan_map[k1].second;
570 : }
571 : else
572 : {
573 830 : _chan_pairs_sf[k].first = _gap_to_chan_map[k1].first;
574 : }
575 :
576 996 : 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 996 : _chan_pairs_sf[k].second = _gap_to_chan_map[k2].first;
583 : }
584 : }
585 : }
586 :
587 : // set the _gij_map
588 4144 : for (unsigned int iz = 0; iz < _n_cells + 1; iz++)
589 : {
590 352008 : for (unsigned int i_gap = 0; i_gap < _n_gaps; i_gap++)
591 : {
592 348030 : if (_gap_type[i_gap] == EChannelType::CENTER)
593 : {
594 266436 : _gij_map[iz].push_back(_pitch - _pin_diameter);
595 : }
596 81594 : else if (_gap_type[i_gap] == EChannelType::EDGE || _gap_type[i_gap] == EChannelType::CORNER)
597 : {
598 81594 : _gij_map[iz].push_back(_duct_to_pin_gap);
599 : }
600 : }
601 : }
602 :
603 11158 : for (unsigned int i = 0; i < _n_channels; i++)
604 : {
605 10992 : if (_subch_type[i] == EChannelType::CENTER || _subch_type[i] == EChannelType::EDGE)
606 : {
607 39984 : for (unsigned int k = 0; k < 3; k++)
608 : {
609 5633028 : for (unsigned int j = 0; j < _n_gaps; j++)
610 : {
611 5603040 : if (_chan_to_gap_map[i][k] == j && i == _gap_to_chan_map[j].first)
612 : {
613 15824 : if (i > _gap_to_chan_map[j].second)
614 : {
615 0 : _sign_id_crossflow_map[i][k] = negative_flow;
616 : }
617 : else
618 : {
619 15824 : _sign_id_crossflow_map[i][k] = positive_flow;
620 : }
621 : }
622 5587216 : else if (_chan_to_gap_map[i][k] == j && i == _gap_to_chan_map[j].second)
623 : {
624 14164 : if (i > _gap_to_chan_map[j].first)
625 : {
626 14164 : _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 996 : else if (_subch_type[i] == EChannelType::CORNER)
637 : {
638 2988 : for (unsigned int k = 0; k < 2; k++)
639 : {
640 193872 : for (unsigned int j = 0; j < _n_gaps; j++)
641 : {
642 191880 : if (_chan_to_gap_map[i][k] == j && i == _gap_to_chan_map[j].first)
643 : {
644 166 : if (i > _gap_to_chan_map[j].second)
645 : {
646 0 : _sign_id_crossflow_map[i][k] = negative_flow;
647 : }
648 : else
649 : {
650 166 : _sign_id_crossflow_map[i][k] = positive_flow;
651 : }
652 : }
653 191714 : else if (_chan_to_gap_map[i][k] == j && i == _gap_to_chan_map[j].second)
654 : {
655 1826 : if (i > _gap_to_chan_map[j].first)
656 : {
657 1826 : _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 11158 : for (unsigned int i = 0; i < _n_channels; i++)
671 : {
672 10992 : if (_subch_type[i] == EChannelType::CENTER)
673 : {
674 7494 : _subchannel_position[i][0] =
675 7494 : (_pin_position[_chan_to_pin_map[i][0]](0) + _pin_position[_chan_to_pin_map[i][1]](0) +
676 7494 : _pin_position[_chan_to_pin_map[i][2]](0)) /
677 : 3.0;
678 7494 : _subchannel_position[i][1] =
679 7494 : (_pin_position[_chan_to_pin_map[i][0]](1) + _pin_position[_chan_to_pin_map[i][1]](1) +
680 7494 : _pin_position[_chan_to_pin_map[i][2]](1)) /
681 : 3.0;
682 : }
683 3498 : else if (_subch_type[i] == EChannelType::EDGE)
684 : {
685 233226 : for (unsigned int j = 0; j < _n_channels; j++)
686 : {
687 230724 : if (_subch_type[j] == EChannelType::CENTER &&
688 170748 : ((_chan_to_pin_map[i][0] == _chan_to_pin_map[j][0] &&
689 2502 : _chan_to_pin_map[i][1] == _chan_to_pin_map[j][1]) ||
690 168246 : (_chan_to_pin_map[i][0] == _chan_to_pin_map[j][1] &&
691 2502 : _chan_to_pin_map[i][1] == _chan_to_pin_map[j][0])))
692 : {
693 2502 : x0 = _pin_position[_chan_to_pin_map[j][2]](0);
694 2502 : y0 = _pin_position[_chan_to_pin_map[j][2]](1);
695 : }
696 228222 : else if (_subch_type[j] == EChannelType::CENTER &&
697 168246 : ((_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 168246 : (_chan_to_pin_map[i][0] == _chan_to_pin_map[j][2] &&
700 1506 : _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 228222 : else if (_subch_type[j] == EChannelType::CENTER &&
706 168246 : ((_chan_to_pin_map[i][0] == _chan_to_pin_map[j][1] &&
707 2502 : _chan_to_pin_map[i][1] == _chan_to_pin_map[j][2]) ||
708 168246 : (_chan_to_pin_map[i][0] == _chan_to_pin_map[j][2] &&
709 1506 : _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 230724 : x1 = 0.5 *
715 230724 : (_pin_position[_chan_to_pin_map[i][0]](0) + _pin_position[_chan_to_pin_map[i][1]](0));
716 230724 : y1 = 0.5 *
717 230724 : (_pin_position[_chan_to_pin_map[i][0]](1) + _pin_position[_chan_to_pin_map[i][1]](1));
718 230724 : a1 = _pin_diameter / 2.0 + _duct_to_pin_gap / 2.0;
719 230724 : a2 = std::sqrt((x1 - x0) * (x1 - x0) + (y1 - y0) * (y1 - y0)) + a1;
720 230724 : _subchannel_position[i][0] = (a2 * x1 - a1 * x0) / (a2 - a1);
721 230724 : _subchannel_position[i][1] = (a2 * y1 - a1 * y0) / (a2 - a1);
722 : } // j
723 : }
724 996 : else if (_subch_type[i] == EChannelType::CORNER)
725 : {
726 996 : x0 = _pin_position[0](0);
727 996 : y0 = _pin_position[0](1);
728 996 : x1 = _pin_position[_chan_to_pin_map[i][0]](0);
729 996 : y1 = _pin_position[_chan_to_pin_map[i][0]](1);
730 996 : a1 = _pin_diameter / 2.0 + _duct_to_pin_gap / 2.0;
731 996 : a2 = std::sqrt((x1 - x0) * (x1 - x0) + (y1 - y0) * (y1 - y0)) + a1;
732 996 : _subchannel_position[i][0] = (a2 * x1 - a1 * x0) / (a2 - a1);
733 996 : _subchannel_position[i][1] = (a2 * y1 - a1 * y0) / (a2 - a1);
734 : }
735 : }
736 :
737 : /// Special case _n_rings == 1
738 166 : 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 11158 : for (auto & gap : _chan_to_gap_map)
751 : {
752 : gap.shrink_to_fit();
753 : }
754 166 : }
755 :
756 : std::unique_ptr<MeshBase>
757 166 : SCMTriSubChannelMeshGenerator::generate()
758 : {
759 166 : auto mesh_base = buildMeshBaseObject();
760 :
761 : BoundaryInfo & boundary_info = mesh_base->get_boundary_info();
762 166 : mesh_base->set_spatial_dimension(3);
763 166 : mesh_base->reserve_elem(_n_cells * _n_channels);
764 166 : mesh_base->reserve_nodes((_n_cells + 1) * _n_channels);
765 166 : _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 11158 : for (unsigned int i = 0; i < _n_channels; i++)
772 : {
773 10992 : _nodes[i].reserve(_n_cells + 1);
774 250968 : for (unsigned int iz = 0; iz < _n_cells + 1; iz++)
775 : {
776 239976 : _nodes[i].push_back(mesh_base->add_point(
777 479952 : 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 11158 : for (unsigned int i = 0; i < _n_channels; i++)
785 : {
786 239976 : for (unsigned int iz = 0; iz < _n_cells; iz++)
787 : {
788 228984 : Elem * elem = new Edge2;
789 228984 : elem->set_id(elem_id++);
790 228984 : elem = mesh_base->add_elem(elem);
791 228984 : const int indx1 = (_n_cells + 1) * i + iz;
792 228984 : const int indx2 = (_n_cells + 1) * i + (iz + 1);
793 228984 : elem->set_node(0, mesh_base->node_ptr(indx1));
794 228984 : elem->set_node(1, mesh_base->node_ptr(indx2));
795 :
796 228984 : if (iz == 0)
797 10992 : boundary_info.add_side(elem, 0, 0);
798 228984 : if (iz == _n_cells - 1)
799 10992 : boundary_info.add_side(elem, 1, 1);
800 : }
801 : }
802 166 : boundary_info.sideset_name(0) = "inlet";
803 166 : boundary_info.sideset_name(1) = "outlet";
804 166 : boundary_info.nodeset_name(0) = "inlet";
805 166 : boundary_info.nodeset_name(1) = "outlet";
806 :
807 : // Naming the block
808 166 : mesh_base->subdomain_name(_block_id) = name();
809 :
810 166 : mesh_base->prepare_for_use();
811 :
812 : // move the meta data into TriSubChannelMesh
813 166 : auto & sch_mesh = static_cast<TriSubChannelMesh &>(*_mesh);
814 166 : sch_mesh._unheated_length_entry = _unheated_length_entry;
815 166 : sch_mesh._heated_length = _heated_length;
816 166 : sch_mesh._unheated_length_exit = _unheated_length_exit;
817 166 : sch_mesh._z_grid = _z_grid;
818 166 : sch_mesh._k_grid = _k_grid;
819 166 : sch_mesh._spacer_z = _spacer_z;
820 166 : sch_mesh._spacer_k = _spacer_k;
821 166 : sch_mesh._z_blockage = _z_blockage;
822 166 : sch_mesh._index_blockage = _index_blockage;
823 166 : sch_mesh._reduction_blockage = _reduction_blockage;
824 166 : sch_mesh._kij = _kij;
825 166 : sch_mesh._pitch = _pitch;
826 166 : sch_mesh._pin_diameter = _pin_diameter;
827 166 : sch_mesh._n_cells = _n_cells;
828 166 : sch_mesh._n_rings = _n_rings;
829 166 : sch_mesh._n_channels = _n_channels;
830 166 : sch_mesh._flat_to_flat = _flat_to_flat;
831 166 : sch_mesh._dwire = _dwire;
832 166 : sch_mesh._hwire = _hwire;
833 166 : sch_mesh._duct_to_pin_gap = _duct_to_pin_gap;
834 166 : sch_mesh._nodes = _nodes;
835 166 : sch_mesh._gap_to_chan_map = _gap_to_chan_map;
836 166 : sch_mesh._gap_to_pin_map = _gap_to_pin_map;
837 166 : sch_mesh._chan_to_gap_map = _chan_to_gap_map;
838 166 : sch_mesh._sign_id_crossflow_map = _sign_id_crossflow_map;
839 166 : sch_mesh._gij_map = _gij_map;
840 166 : sch_mesh._subchannel_position = _subchannel_position;
841 166 : sch_mesh._pin_position = _pin_position;
842 166 : sch_mesh._pins_in_rings = _pins_in_rings;
843 166 : sch_mesh._chan_to_pin_map = _chan_to_pin_map;
844 166 : sch_mesh._npins = _npins;
845 166 : sch_mesh._n_gaps = _n_gaps;
846 166 : sch_mesh._subch_type = _subch_type;
847 166 : sch_mesh._gap_type = _gap_type;
848 166 : sch_mesh._gap_pairs_sf = _gap_pairs_sf;
849 166 : sch_mesh._chan_pairs_sf = _chan_pairs_sf;
850 166 : sch_mesh._pin_to_chan_map = _pin_to_chan_map;
851 :
852 166 : return mesh_base;
853 0 : }
|
