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