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 "BreakMeshByBlockGenerator.h"
11 : #include "CastUniquePointer.h"
12 : #include "MooseMeshUtils.h"
13 :
14 : #include "libmesh/distributed_mesh.h"
15 : #include "libmesh/elem.h"
16 : #include "libmesh/partitioner.h"
17 : #include <typeinfo>
18 :
19 : registerMooseObject("MooseApp", BreakMeshByBlockGenerator);
20 :
21 : InputParameters
22 14737 : BreakMeshByBlockGenerator::validParams()
23 : {
24 14737 : InputParameters params = BreakMeshByBlockGeneratorBase::validParams();
25 14737 : params.addRequiredParam<MeshGeneratorName>("input", "The mesh we want to modify");
26 14737 : params.addClassDescription(
27 : "Break the mesh at interfaces between blocks. New nodes will be generated so elements on "
28 : "each side of the break are no longer connected. At the moment, this only works on a "
29 : "REPLICATED mesh");
30 14737 : params.addParam<std::vector<SubdomainName>>(
31 : "surrounding_blocks",
32 : "The list of subdomain names surrounding which interfaces will be generated.");
33 14737 : params.addParam<std::vector<std::vector<SubdomainName>>>(
34 : "block_pairs", "The list of subdomain pairs between which interfaces will be generated.");
35 44211 : params.addParam<bool>("add_transition_interface",
36 29474 : false,
37 : "If true and block is not empty, a special boundary named "
38 : "interface_transition is generate between listed blocks and other blocks.");
39 44211 : params.addParam<bool>(
40 29474 : "split_transition_interface", false, "Whether to split the transition interface by blocks.");
41 44211 : params.addParam<bool>("add_interface_on_two_sides",
42 29474 : false,
43 : "Whether to add an additional interface boundary at the other side.");
44 14737 : params.addParam<BoundaryName>(
45 : "interface_transition_name",
46 : "interface_transition",
47 : "the name of the interface transition boundary created when blocks are provided");
48 14737 : return params;
49 0 : }
50 :
51 236 : BreakMeshByBlockGenerator::BreakMeshByBlockGenerator(const InputParameters & parameters)
52 : : BreakMeshByBlockGeneratorBase(parameters),
53 236 : _input(getMesh("input")),
54 236 : _block_pairs_restricted(parameters.isParamSetByUser("block_pairs")),
55 236 : _surrounding_blocks_restricted(parameters.isParamSetByUser("surrounding_blocks")),
56 236 : _add_transition_interface(getParam<bool>("add_transition_interface")),
57 236 : _split_transition_interface(getParam<bool>("split_transition_interface")),
58 236 : _interface_transition_name(getParam<BoundaryName>("interface_transition_name")),
59 708 : _add_interface_on_two_sides(getParam<bool>("add_interface_on_two_sides"))
60 : {
61 236 : if (_block_pairs_restricted && _surrounding_blocks_restricted)
62 0 : paramError("block_pairs_restricted",
63 : "BreakMeshByBlockGenerator: 'surrounding_blocks' and 'block_pairs' can not be used "
64 : "at the same time.");
65 :
66 236 : if (!_add_transition_interface && _split_transition_interface)
67 0 : paramError("split_transition_interface",
68 : "BreakMeshByBlockGenerator cannot split the transition interface because "
69 : "add_transition_interface is false");
70 :
71 236 : if (_add_transition_interface && _block_pairs_restricted)
72 0 : paramError("add_transition_interface",
73 : "BreakMeshByBlockGenerator cannot split the transition interface when 'block_pairs' "
74 : "are specified.");
75 236 : }
76 :
77 : std::unique_ptr<MeshBase>
78 236 : BreakMeshByBlockGenerator::generate()
79 : {
80 236 : std::unique_ptr<MeshBase> mesh = std::move(_input);
81 236 : if (!mesh->is_replicated())
82 0 : mooseError("BreakMeshByBlockGenerator is not implemented for distributed meshes");
83 :
84 : // Try to trick the rest of the world into thinking we're prepared
85 236 : mesh->prepare_for_use();
86 :
87 236 : BoundaryInfo & boundary_info = mesh->get_boundary_info();
88 :
89 : // Handle block restrictions
90 236 : if (_block_pairs_restricted)
91 : {
92 60 : for (const auto & block_name_pair :
93 200 : getParam<std::vector<std::vector<SubdomainName>>>("block_pairs"))
94 : {
95 84 : if (block_name_pair.size() != 2)
96 0 : paramError("block_pairs",
97 : "Each row of 'block_pairs' must have a size of two (block names).");
98 :
99 : // check that the blocks exist in the mesh
100 248 : for (const auto & name : block_name_pair)
101 168 : if (!MooseMeshUtils::hasSubdomainName(*mesh, name))
102 4 : paramError("block_pairs", "The block '", name, "' was not found in the mesh");
103 :
104 80 : const auto block_pair = MooseMeshUtils::getSubdomainIDs(*mesh, block_name_pair);
105 160 : std::pair<SubdomainID, SubdomainID> pair = std::make_pair(
106 160 : std::min(block_pair[0], block_pair[1]), std::max(block_pair[0], block_pair[1]));
107 :
108 80 : _block_pairs.insert(pair);
109 80 : std::copy(block_pair.begin(), block_pair.end(), std::inserter(_block_set, _block_set.end()));
110 80 : }
111 : }
112 176 : else if (_surrounding_blocks_restricted)
113 : {
114 : // check that the blocks exist in the mesh
115 1284 : for (const auto & name : getParam<std::vector<SubdomainName>>("surrounding_blocks"))
116 1244 : if (!MooseMeshUtils::hasSubdomainName(*mesh, name))
117 4 : paramError("surrounding_blocks", "The block '", name, "' was not found in the mesh");
118 :
119 : const auto surrounding_block_ids = MooseMeshUtils::getSubdomainIDs(
120 40 : *mesh, getParam<std::vector<SubdomainName>>("surrounding_blocks"));
121 80 : std::copy(surrounding_block_ids.begin(),
122 : surrounding_block_ids.end(),
123 40 : std::inserter(_block_set, _block_set.end()));
124 40 : }
125 :
126 : // check that a boundary named _interface_transition_name does not already exist in the mesh
127 252 : if ((_block_pairs_restricted || _surrounding_blocks_restricted) && _add_transition_interface &&
128 24 : boundary_info.get_id_by_name(_interface_transition_name) != Moose::INVALID_BOUNDARY_ID)
129 0 : paramError("interface_transition_name",
130 : "BreakMeshByBlockGenerator the specified interface transition boundary name "
131 : "already exits");
132 :
133 : // initialize the node to element map
134 228 : std::map<dof_id_type, std::vector<dof_id_type>> node_to_elem_map;
135 37044 : for (const auto & elem : mesh->active_element_ptr_range())
136 106400 : for (unsigned int n = 0; n < elem->n_nodes(); n++)
137 88220 : node_to_elem_map[elem->node_id(n)].push_back(elem->id());
138 :
139 16148 : for (auto node_it = node_to_elem_map.begin(); node_it != node_to_elem_map.end(); ++node_it)
140 : {
141 15920 : const dof_id_type current_node_id = node_it->first;
142 15920 : const Node * current_node = mesh->node_ptr(current_node_id);
143 :
144 15920 : if (current_node != nullptr)
145 : {
146 : // find node multiplicity
147 15920 : std::set<subdomain_id_type> connected_blocks;
148 103912 : for (auto elem_id = node_it->second.begin(); elem_id != node_it->second.end(); elem_id++)
149 : {
150 87992 : const Elem * current_elem = mesh->elem_ptr(*elem_id);
151 87992 : subdomain_id_type block_id = blockRestrictedElementSubdomainID(current_elem);
152 87992 : if (!_block_pairs_restricted)
153 81976 : connected_blocks.insert(block_id);
154 : else
155 : {
156 6016 : if (block_id != Elem::invalid_subdomain_id)
157 5760 : connected_blocks.insert(block_id);
158 : }
159 : }
160 :
161 15920 : unsigned int node_multiplicity = connected_blocks.size();
162 :
163 : // check if current_node need to be duplicated
164 15920 : if (node_multiplicity > 1)
165 : {
166 : // retrieve connected elements from the map
167 7676 : const std::vector<dof_id_type> & connected_elems = node_it->second;
168 :
169 : // find reference_subdomain_id (e.g. the subdomain with lower id or the
170 : // Elem::invalid_subdomain_id)
171 7676 : auto subdomain_it = connected_blocks.begin();
172 : subdomain_id_type reference_subdomain_id =
173 7676 : connected_blocks.find(Elem::invalid_subdomain_id) != connected_blocks.end()
174 13192 : ? Elem::invalid_subdomain_id
175 5516 : : *subdomain_it;
176 :
177 : // For the block_pairs option, the node that is only shared by specific block pairs will be
178 : // newly created.
179 7676 : bool should_create_new_node = true;
180 7676 : if (_block_pairs_restricted)
181 : {
182 456 : auto elems = node_to_elem_map[current_node->id()];
183 456 : should_create_new_node = false;
184 456 : std::set<subdomain_id_type> sets_blocks_for_this_node;
185 2216 : for (auto elem_id = elems.begin(); elem_id != elems.end(); elem_id++)
186 1760 : sets_blocks_for_this_node.insert(
187 1760 : blockRestrictedElementSubdomainID(mesh->elem_ptr(*elem_id)));
188 456 : if (sets_blocks_for_this_node.size() == 2)
189 : {
190 408 : auto setIt = sets_blocks_for_this_node.begin();
191 408 : if (findBlockPairs(*setIt, *std::next(setIt, 1)))
192 352 : should_create_new_node = true;
193 : }
194 456 : }
195 :
196 : // multiplicity counter to keep track of how many nodes we added
197 7676 : unsigned int multiplicity_counter = node_multiplicity;
198 60020 : for (auto elem_id : connected_elems)
199 : {
200 : // all the duplicate nodes are added and assigned
201 52352 : if (multiplicity_counter == 0)
202 8 : break;
203 :
204 52344 : Elem * current_elem = mesh->elem_ptr(elem_id);
205 52344 : subdomain_id_type block_id = blockRestrictedElementSubdomainID(current_elem);
206 :
207 75040 : if ((blockRestrictedElementSubdomainID(current_elem) != reference_subdomain_id) ||
208 22696 : (_block_pairs_restricted & findBlockPairs(reference_subdomain_id, block_id)))
209 : {
210 : // assign the newly added node to current_elem
211 29648 : Node * new_node = nullptr;
212 :
213 29648 : std::vector<boundary_id_type> node_boundary_ids;
214 :
215 141120 : for (unsigned int node_id = 0; node_id < current_elem->n_nodes(); ++node_id)
216 126796 : if (current_elem->node_id(node_id) ==
217 126796 : current_node->id()) // if current node == node on element
218 : {
219 15324 : if (should_create_new_node)
220 : {
221 : // add new node
222 15164 : new_node = Node::build(*current_node, mesh->n_nodes()).release();
223 :
224 : // We're duplicating nodes so that each subdomain elem has its own copy, so it
225 : // seems natural to assign this new node the same proc id as corresponding
226 : // subdomain elem
227 15164 : new_node->processor_id() = current_elem->processor_id();
228 15164 : mesh->add_node(new_node);
229 15164 : current_elem->set_node(node_id, new_node);
230 : // Add boundary info to the new node
231 15164 : boundary_info.boundary_ids(current_node, node_boundary_ids);
232 15164 : boundary_info.add_node(new_node, node_boundary_ids);
233 : }
234 :
235 15324 : multiplicity_counter--; // node created, update multiplicity counter
236 :
237 15324 : break; // ones the proper node has been fixed in one element we can break the
238 : // loop
239 : }
240 :
241 29648 : if (should_create_new_node)
242 : {
243 419248 : for (auto connected_elem_id : connected_elems)
244 : {
245 389760 : Elem * connected_elem = mesh->elem_ptr(connected_elem_id);
246 :
247 : // Assign the newly added node to other connected elements with the same
248 : // block_id
249 389760 : if (connected_elem->subdomain_id() == current_elem->subdomain_id() &&
250 : connected_elem != current_elem)
251 : {
252 633832 : for (unsigned int node_id = 0; node_id < connected_elem->n_nodes(); ++node_id)
253 516220 : if (connected_elem->node_id(node_id) ==
254 516220 : current_node->id()) // if current node == node on element
255 : {
256 14324 : connected_elem->set_node(node_id, new_node);
257 14324 : break;
258 : }
259 : }
260 : }
261 : }
262 29648 : }
263 : }
264 :
265 : // create blocks pair and assign element side to new interface boundary map
266 60036 : for (auto elem_id : connected_elems)
267 : {
268 741648 : for (auto connected_elem_id : connected_elems)
269 : {
270 689288 : Elem * current_elem = mesh->elem_ptr(elem_id);
271 689288 : Elem * connected_elem = mesh->elem_ptr(connected_elem_id);
272 689288 : subdomain_id_type curr_elem_subid = blockRestrictedElementSubdomainID(current_elem);
273 : subdomain_id_type connected_elem_subid =
274 689288 : blockRestrictedElementSubdomainID(connected_elem);
275 :
276 689288 : if (current_elem != connected_elem && curr_elem_subid < connected_elem_subid)
277 : {
278 193536 : if (current_elem->has_neighbor(connected_elem))
279 : {
280 30032 : dof_id_type elem_id = current_elem->id();
281 30032 : dof_id_type connected_elem_id = connected_elem->id();
282 30032 : unsigned int side = current_elem->which_neighbor_am_i(connected_elem);
283 : unsigned int connected_elem_side =
284 30032 : connected_elem->which_neighbor_am_i(current_elem);
285 :
286 : // in this case we need to play a game to reorder the sides
287 30032 : if (_block_pairs_restricted || _surrounding_blocks_restricted)
288 : {
289 16128 : connected_elem_subid = connected_elem->subdomain_id();
290 16128 : if (curr_elem_subid > connected_elem_subid) // we need to switch the ids
291 : {
292 2896 : connected_elem_subid = current_elem->subdomain_id();
293 2896 : curr_elem_subid = connected_elem->subdomain_id();
294 :
295 2896 : elem_id = connected_elem->id();
296 2896 : side = connected_elem->which_neighbor_am_i(current_elem);
297 :
298 2896 : connected_elem_id = current_elem->id();
299 2896 : connected_elem_side = current_elem->which_neighbor_am_i(connected_elem);
300 : }
301 : }
302 :
303 : std::pair<subdomain_id_type, subdomain_id_type> blocks_pair =
304 30032 : std::make_pair(curr_elem_subid, connected_elem_subid);
305 :
306 : std::pair<subdomain_id_type, subdomain_id_type> blocks_pair2 =
307 30032 : std::make_pair(connected_elem_subid, curr_elem_subid);
308 :
309 30032 : if (_block_pairs_restricted)
310 : {
311 928 : if (findBlockPairs(blockRestrictedElementSubdomainID(current_elem),
312 928 : blockRestrictedElementSubdomainID(connected_elem)))
313 : {
314 768 : _new_boundary_sides_map[blocks_pair].insert(std::make_pair(elem_id, side));
315 768 : if (_add_interface_on_two_sides)
316 1216 : _new_boundary_sides_map[blocks_pair2].insert(
317 1216 : std::make_pair(connected_elem_id, connected_elem_side));
318 : }
319 : }
320 : else
321 : {
322 29104 : _new_boundary_sides_map[blocks_pair].insert(std::make_pair(elem_id, side));
323 29104 : if (_add_interface_on_two_sides)
324 0 : _new_boundary_sides_map[blocks_pair2].insert(
325 0 : std::make_pair(connected_elem_id, connected_elem_side));
326 : }
327 : }
328 : }
329 : }
330 : }
331 :
332 : } // end multiplicity check
333 15920 : } // end loop over nodes
334 : } // end nodeptr check
335 :
336 228 : addInterfaceBoundary(*mesh);
337 228 : Partitioner::set_node_processor_ids(*mesh);
338 456 : return dynamic_pointer_cast<MeshBase>(mesh);
339 228 : }
340 :
341 : void
342 228 : BreakMeshByBlockGenerator::addInterfaceBoundary(MeshBase & mesh)
343 : {
344 228 : BoundaryInfo & boundary_info = mesh.get_boundary_info();
345 :
346 : boundary_id_type boundary_id;
347 228 : boundary_id_type boundary_id_interface = Moose::INVALID_BOUNDARY_ID;
348 228 : boundary_id_type boundary_id_interface_transition = Moose::INVALID_BOUNDARY_ID;
349 :
350 228 : BoundaryName boundary_name;
351 :
352 : // loop over boundary sides
353 5004 : for (auto & boundary_side_map : _new_boundary_sides_map)
354 : {
355 4776 : if (!(_block_pairs_restricted || _surrounding_blocks_restricted) ||
356 3912 : ((_block_pairs_restricted || _surrounding_blocks_restricted) && !_add_transition_interface))
357 : {
358 : // find the appropriate boundary name and id
359 : // given primary and secondary block ID
360 2496 : if (_split_interface)
361 1320 : findBoundaryNameAndInd(mesh,
362 1320 : boundary_side_map.first.first,
363 1320 : boundary_side_map.first.second,
364 : boundary_name,
365 : boundary_id,
366 : boundary_info);
367 : else
368 : {
369 1176 : boundary_name = _interface_name;
370 2352 : boundary_id_interface = boundary_id_interface == Moose::INVALID_BOUNDARY_ID
371 1176 : ? findFreeBoundaryId(mesh)
372 : : boundary_id_interface;
373 1176 : boundary_id = boundary_id_interface;
374 1176 : boundary_info.sideset_name(boundary_id_interface) = boundary_name;
375 : }
376 : }
377 : else // block resticted with transition boundary
378 : {
379 :
380 : const bool interior_boundary =
381 4128 : _block_set.find(boundary_side_map.first.first) != _block_set.end() &&
382 1848 : _block_set.find(boundary_side_map.first.second) != _block_set.end();
383 :
384 2280 : if ((_split_interface && interior_boundary) ||
385 1320 : (_split_transition_interface && !interior_boundary))
386 : {
387 1520 : findBoundaryNameAndInd(mesh,
388 1520 : boundary_side_map.first.first,
389 1520 : boundary_side_map.first.second,
390 : boundary_name,
391 : boundary_id,
392 : boundary_info);
393 : }
394 760 : else if (interior_boundary)
395 : {
396 480 : boundary_name = _interface_name;
397 960 : boundary_id_interface = boundary_id_interface == Moose::INVALID_BOUNDARY_ID
398 480 : ? findFreeBoundaryId(mesh)
399 : : boundary_id_interface;
400 :
401 480 : boundary_id = boundary_id_interface;
402 480 : boundary_info.sideset_name(boundary_id_interface) = boundary_name;
403 : }
404 : else
405 : {
406 280 : boundary_id_interface_transition =
407 280 : boundary_id_interface_transition == Moose::INVALID_BOUNDARY_ID
408 280 : ? findFreeBoundaryId(mesh)
409 : : boundary_id_interface_transition;
410 280 : boundary_id = boundary_id_interface_transition;
411 280 : boundary_info.sideset_name(boundary_id) = _interface_transition_name;
412 : }
413 : }
414 : // loop over all the side belonging to each pair and add it to the proper interface
415 13392 : for (auto & element_side : boundary_side_map.second)
416 8616 : boundary_info.add_side(element_side.first, element_side.second, boundary_id);
417 : }
418 228 : }
419 :
420 : subdomain_id_type
421 1574872 : BreakMeshByBlockGenerator::blockRestrictedElementSubdomainID(const Elem * elem)
422 : {
423 1574872 : subdomain_id_type elem_subdomain_id = elem->subdomain_id();
424 1844632 : if ((_block_pairs_restricted || _surrounding_blocks_restricted) &&
425 1844632 : (_block_set.find(elem_subdomain_id) == _block_set.end()))
426 88672 : elem_subdomain_id = Elem::invalid_subdomain_id;
427 :
428 1574872 : return elem_subdomain_id;
429 : }
430 :
431 : bool
432 24032 : BreakMeshByBlockGenerator::findBlockPairs(subdomain_id_type block_one, subdomain_id_type block_two)
433 : {
434 25592 : for (auto block_pair : _block_pairs)
435 2680 : if ((block_pair.first == block_one && block_pair.second == block_two) ||
436 1560 : (block_pair.first == block_two && block_pair.second == block_one))
437 1120 : return true;
438 22912 : return false;
439 : }
|
