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