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 "SCMQuadDuctMeshGenerator.h"
11 : #include "QuadSubChannelMesh.h"
12 : #include <cmath>
13 : #include "libmesh/unstructured_mesh.h"
14 : #include "libmesh/face_quad4.h"
15 :
16 : registerMooseObject("SubChannelApp", SCMQuadDuctMeshGenerator);
17 :
18 : InputParameters
19 10 : SCMQuadDuctMeshGenerator::validParams()
20 : {
21 10 : InputParameters params = MeshGenerator::validParams();
22 10 : params.addClassDescription("Creates a mesh of 2D duct cells around a square-lattice subassembly");
23 20 : params.addRequiredParam<MeshGeneratorName>("input", "The corresponding subchannel mesh");
24 20 : params.addParam<unsigned int>("block_id", 2, "Subdomain id for the duct mesh cells");
25 20 : params.addRequiredParam<unsigned int>("n_cells", "The number of cells in the axial direction");
26 20 : params.addParam<Real>("unheated_length_entry", 0.0, "Unheated length at entry [m]");
27 20 : params.addRequiredParam<Real>("heated_length", "Heated length [m]");
28 20 : params.addParam<Real>("unheated_length_exit", 0.0, "Unheated length at exit [m]");
29 20 : params.addRangeCheckedParam<Real>("pitch", "pitch > 0", "Lattice pitch (must be positive)");
30 20 : params.addRangeCheckedParam<unsigned int>(
31 : "nx",
32 : "nx > 1",
33 : "Number of channels in the x direction for the subchannel assembly. Must be more than 1 to "
34 : "built a duct[-]");
35 20 : params.addRangeCheckedParam<unsigned int>(
36 : "ny",
37 : "ny > 1",
38 : "Number of channels in the y direction for the subchannel assembly. Must be more than 1 to "
39 : "built a duct[-]");
40 20 : params.addRequiredParam<Real>("side_gap",
41 : "Gap between duct wall and outer pin lattice: distance(edge pin "
42 : "center, duct wall) = pitch/2 + side_gap [m]");
43 10 : return params;
44 0 : }
45 :
46 5 : SCMQuadDuctMeshGenerator::SCMQuadDuctMeshGenerator(const InputParameters & params)
47 : : MeshGenerator(params),
48 5 : _input(getMesh("input")),
49 10 : _n_cells(getParam<unsigned int>("n_cells")),
50 10 : _unheated_length_entry(getParam<Real>("unheated_length_entry")),
51 10 : _heated_length(getParam<Real>("heated_length")),
52 10 : _unheated_length_exit(getParam<Real>("unheated_length_exit")),
53 10 : _block_id(getParam<unsigned int>("block_id")),
54 10 : _pitch(getParam<Real>("pitch")),
55 10 : _nx(getParam<unsigned int>("nx")),
56 10 : _ny(getParam<unsigned int>("ny")),
57 15 : _side_gap(getParam<Real>("side_gap"))
58 : {
59 5 : SubChannelMesh::generateZGrid(
60 5 : _unheated_length_entry, _heated_length, _unheated_length_exit, _n_cells, _z_grid);
61 5 : }
62 :
63 : std::unique_ptr<MeshBase>
64 5 : SCMQuadDuctMeshGenerator::generate()
65 : {
66 5 : std::unique_ptr<MeshBase> mesh_base = std::move(_input);
67 5 : mesh_base->set_mesh_dimension(3);
68 :
69 : std::vector<Point> cross_sec;
70 5 : ductCrossSec(cross_sec, _nx, _ny, _pitch, _side_gap);
71 : std::vector<Point> points;
72 5 : ductPoints(points, cross_sec, _z_grid);
73 : std::vector<std::vector<size_t>> elem_point_indices;
74 5 : ductElems(elem_point_indices, _z_grid.size(), cross_sec.size());
75 : std::vector<Node *> duct_nodes;
76 5 : buildDuct(mesh_base, duct_nodes, points, elem_point_indices, _block_id);
77 5 : mesh_base->subdomain_name(_block_id) = name();
78 :
79 5 : mesh_base->prepare_for_use();
80 :
81 : // Mirror the Tri variant: provide mapping hooks into the subchannel mesh
82 5 : auto & sch_mesh = static_cast<QuadSubChannelMesh &>(*_mesh);
83 5 : sch_mesh.setChannelToDuctMaps(duct_nodes);
84 :
85 5 : return mesh_base;
86 5 : }
87 :
88 : size_t
89 3570 : SCMQuadDuctMeshGenerator::ductPointIndex(unsigned int points_per_layer,
90 : unsigned int layer,
91 : unsigned int point) const
92 : {
93 3570 : return layer * points_per_layer + point;
94 : }
95 :
96 : void
97 0 : SCMQuadDuctMeshGenerator::ductCorners(std::vector<Point> & corners,
98 : Real half_x,
99 : Real half_y,
100 : const Point & center) const
101 : {
102 0 : corners.resize(4);
103 0 : corners[0] = center + Point(-half_x, -half_y, 0);
104 0 : corners[1] = center + Point(half_x, -half_y, 0);
105 0 : corners[2] = center + Point(half_x, half_y, 0);
106 0 : corners[3] = center + Point(-half_x, half_y, 0);
107 0 : }
108 :
109 : void
110 5 : SCMQuadDuctMeshGenerator::ductCrossSec(std::vector<Point> & cross_sec,
111 : unsigned int nx,
112 : unsigned int ny,
113 : Real pitch,
114 : Real side_gap) const
115 : {
116 5 : cross_sec.clear();
117 :
118 5 : const Real half_x = 0.5 * (nx - 1) * pitch + side_gap;
119 5 : const Real half_y = 0.5 * (ny - 1) * pitch + side_gap;
120 :
121 : // Subchannel-aligned grid extents (no side_gap)
122 : const Real x0 = -0.5 * (nx - 1) * pitch;
123 : const Real y0 = -0.5 * (ny - 1) * pitch;
124 :
125 : // ---- Bottom edge (y = -half_y): nx points, left -> right ----
126 30 : for (unsigned int i = 0; i < nx; ++i)
127 : {
128 : Real x;
129 25 : if (i == 0)
130 5 : x = -half_x;
131 20 : else if (i == nx - 1)
132 5 : x = half_x;
133 : else
134 15 : x = x0 + i * pitch; // aligned with subchannel x grid
135 :
136 25 : cross_sec.emplace_back(x, -half_y, 0.0);
137 : }
138 :
139 : // ---- Right edge (x = +half_x): (ny-2) points, bottom -> top (exclude corners) ----
140 15 : for (unsigned int j = 1; j + 1 < ny; ++j)
141 : {
142 10 : const Real y = y0 + j * pitch; // aligned with subchannel y grid
143 10 : cross_sec.emplace_back(half_x, y, 0.0);
144 : }
145 :
146 : // ---- Top edge (y = +half_y): nx points, right -> left ----
147 30 : for (unsigned int i = 0; i < nx; ++i)
148 : {
149 25 : const unsigned int ii = nx - 1 - i;
150 :
151 : Real x;
152 25 : if (ii == 0)
153 5 : x = -half_x;
154 20 : else if (ii == nx - 1)
155 5 : x = half_x;
156 : else
157 15 : x = x0 + ii * pitch;
158 :
159 25 : cross_sec.emplace_back(x, half_y, 0.0);
160 : }
161 :
162 : // ---- Left edge (x = -half_x): (ny-2) points, top -> bottom (exclude corners) ----
163 15 : for (unsigned int j = 1; j + 1 < ny; ++j)
164 : {
165 10 : const unsigned int jj = ny - 1 - j;
166 10 : const Real y = y0 + jj * pitch; // aligned with subchannel y grid
167 10 : cross_sec.emplace_back(-half_x, y, 0.0);
168 : }
169 :
170 : // Total points: 2*nx + 2*ny - 4
171 5 : }
172 :
173 : void
174 5 : SCMQuadDuctMeshGenerator::ductPoints(std::vector<Point> & points,
175 : const std::vector<Point> & cross_sec,
176 : const std::vector<Real> & z_layers) const
177 : {
178 5 : points.resize(cross_sec.size() * z_layers.size());
179 60 : for (size_t i = 0; i < z_layers.size(); i++)
180 825 : for (size_t j = 0; j < cross_sec.size(); j++)
181 770 : points[ductPointIndex(cross_sec.size(), i, j)] =
182 770 : Point(cross_sec[j](0), cross_sec[j](1), z_layers[i]);
183 5 : }
184 :
185 : void
186 5 : SCMQuadDuctMeshGenerator::ductElems(std::vector<std::vector<size_t>> & elem_point_indices,
187 : unsigned int n_layers,
188 : unsigned int points_per_layer) const
189 : {
190 5 : elem_point_indices.clear();
191 55 : for (unsigned int i = 0; i < n_layers - 1; i++)
192 : {
193 : const unsigned int bottom = i;
194 50 : const unsigned int top = i + 1;
195 750 : for (unsigned int j = 0; j < points_per_layer; j++)
196 : {
197 : const unsigned int left = j;
198 700 : const unsigned int right = (j + 1) % points_per_layer;
199 2100 : elem_point_indices.push_back({ductPointIndex(points_per_layer, bottom, left),
200 700 : ductPointIndex(points_per_layer, bottom, right),
201 700 : ductPointIndex(points_per_layer, top, right),
202 700 : ductPointIndex(points_per_layer, top, left)});
203 : }
204 : }
205 5 : }
206 :
207 : void
208 5 : SCMQuadDuctMeshGenerator::buildDuct(std::unique_ptr<MeshBase> & mesh,
209 : std::vector<Node *> & duct_nodes,
210 : const std::vector<Point> & points,
211 : const std::vector<std::vector<size_t>> & elem_point_indices,
212 : SubdomainID block) const
213 : {
214 : // Create mesh nodes for all duct points and keep a local index -> Node* map
215 5 : duct_nodes.clear();
216 5 : duct_nodes.reserve(points.size());
217 775 : for (const auto & p : points)
218 770 : duct_nodes.push_back(mesh->add_point(p));
219 :
220 : // Create QUAD4 surface elements using libMesh factory style
221 705 : for (const auto & elem_indices : elem_point_indices)
222 : {
223 : mooseAssert(elem_indices.size() == 4,
224 : "Expected 4 node indices per element when building QUAD4 elements.");
225 :
226 700 : auto elem = Elem::build(ElemType::QUAD4);
227 700 : elem->subdomain_id() = block;
228 :
229 : // Set the 4 nodes of the QUAD4
230 3500 : for (unsigned int i = 0; i < 4; ++i)
231 : {
232 2800 : const auto idx = elem_indices[i];
233 : mooseAssert(idx < duct_nodes.size(), "Element node index out of range.");
234 2800 : elem->set_node(i, duct_nodes[idx]);
235 : }
236 :
237 : // Hand ownership to the mesh
238 700 : mesh->add_elem(elem.release());
239 700 : }
240 5 : }
|