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 "BreakMeshByElementGenerator.h"
11 : #include "CastUniquePointer.h"
12 : #include "MooseMeshUtils.h"
13 :
14 : #include "libmesh/partitioner.h"
15 :
16 : registerMooseObject("MooseApp", BreakMeshByElementGenerator);
17 : registerMooseObjectRenamed("MooseApp",
18 : ExplodeMeshGenerator,
19 : "05/18/2024 24:00",
20 : BreakMeshByElementGenerator);
21 :
22 : InputParameters
23 6232 : BreakMeshByElementGenerator::validParams()
24 : {
25 6232 : InputParameters params = MeshGenerator::validParams();
26 12464 : params.addClassDescription("Break all element-element interfaces in the specified subdomains.");
27 24928 : params.addRequiredParam<MeshGeneratorName>("input", "The mesh we want to modify");
28 18696 : params.addParam<std::vector<SubdomainID>>(
29 : "subdomains",
30 12464 : std::vector<SubdomainID>(),
31 : "The list of subdomain IDs to explode. Leave unset to explode all subdomains.");
32 24928 : params.addParam<BoundaryName>(
33 : "interface_name",
34 : "element_boundaries",
35 : "The boundary name containing all broken element-element interfaces.");
36 24928 : params.addRangeCheckedParam<unsigned int>(
37 : "interface_sides",
38 12464 : 1,
39 : "interface_sides<3",
40 : "Whether to add no interface boundary, a 1-sided boundary (facing from lower to higher "
41 : "element id), or a 2-sided boundary");
42 6232 : return params;
43 0 : }
44 :
45 55 : BreakMeshByElementGenerator::BreakMeshByElementGenerator(const InputParameters & parameters)
46 : : MeshGenerator(parameters),
47 55 : _input(getMesh("input")),
48 110 : _subdomains(getParam<std::vector<SubdomainID>>("subdomains")),
49 110 : _interface_name(getParam<BoundaryName>("interface_name")),
50 165 : _interface_sides(getParam<unsigned int>("interface_sides"))
51 : {
52 55 : }
53 :
54 : std::unique_ptr<MeshBase>
55 52 : BreakMeshByElementGenerator::generate()
56 : {
57 52 : std::unique_ptr<MeshBase> mesh = std::move(_input);
58 :
59 52 : if (!mesh->is_prepared())
60 19 : mesh->prepare_for_use();
61 :
62 : // check that the subdomain IDs exist in the mesh
63 106 : for (const auto & id : _subdomains)
64 57 : if (!MooseMeshUtils::hasSubdomainID(*mesh, id))
65 6 : paramError("subdomains", "The block ID '", id, "' was not found in the mesh");
66 :
67 49 : BoundaryInfo & boundary_info = mesh->get_boundary_info();
68 98 : if (_interface_sides &&
69 49 : boundary_info.get_id_by_name(_interface_name) != Moose::INVALID_BOUNDARY_ID)
70 0 : paramError("interface_name", "The specified interface name already exists in the mesh.");
71 :
72 49 : const auto node_to_elem_map = buildSubdomainRestrictedNodeToElemMap(mesh, _subdomains);
73 :
74 49 : duplicateNodes(mesh, node_to_elem_map);
75 :
76 49 : createInterface(*mesh, node_to_elem_map);
77 :
78 49 : Partitioner::set_node_processor_ids(*mesh);
79 :
80 : // We need to update the global_boundary_ids, and this is faster
81 : // than a full prepare_for_use()
82 49 : boundary_info.regenerate_id_sets();
83 :
84 98 : return dynamic_pointer_cast<MeshBase>(mesh);
85 49 : }
86 :
87 : BreakMeshByElementGenerator::NodeToElemMapType
88 49 : BreakMeshByElementGenerator::buildSubdomainRestrictedNodeToElemMap(
89 : std::unique_ptr<MeshBase> & mesh, const std::vector<SubdomainID> & subdomains) const
90 : {
91 49 : NodeToElemMapType node_to_elem_map;
92 54167 : for (const auto & elem : mesh->active_element_ptr_range())
93 : {
94 : // Skip if subdomains are specified and the element is not in them
95 28422 : if (!subdomains.empty() &&
96 28422 : std::find(subdomains.begin(), subdomains.end(), elem->subdomain_id()) == subdomains.end())
97 192 : continue;
98 :
99 26867 : std::set<const Elem *> neighbors;
100 26867 : elem->find_point_neighbors(neighbors);
101 :
102 134847 : for (auto n : make_range(elem->n_nodes()))
103 : {
104 : // if ANY neighboring element that contains this node is not in specified subdomains,
105 : // don't add this node to the map, i.e. don't split this node.
106 107980 : bool should_duplicate = true;
107 107980 : if (!subdomains.empty())
108 43309 : for (auto neighbor : neighbors)
109 55713 : if (neighbor->contains_point(elem->node_ref(n)) &&
110 17280 : std::find(subdomains.begin(), subdomains.end(), neighbor->subdomain_id()) ==
111 55713 : subdomains.end())
112 : {
113 320 : should_duplicate = false;
114 320 : break;
115 : }
116 :
117 107980 : if (should_duplicate)
118 107660 : node_to_elem_map[elem->node_id(n)].insert(elem->id());
119 : }
120 26916 : }
121 :
122 49 : return node_to_elem_map;
123 0 : }
124 :
125 : void
126 49 : BreakMeshByElementGenerator::duplicateNodes(std::unique_ptr<MeshBase> & mesh,
127 : const NodeToElemMapType & node_to_elem_map) const
128 : {
129 8740 : for (const auto & [node_id, connected_elem_ids] : node_to_elem_map)
130 116351 : for (auto & connected_elem_id : connected_elem_ids)
131 107660 : if (connected_elem_id != *connected_elem_ids.begin())
132 98969 : duplicateNode(mesh, mesh->elem_ptr(connected_elem_id), mesh->node_ptr(node_id));
133 49 : }
134 :
135 : void
136 98969 : BreakMeshByElementGenerator::duplicateNode(std::unique_ptr<MeshBase> & mesh,
137 : Elem * elem,
138 : const Node * node) const
139 : {
140 98969 : std::unique_ptr<Node> new_node = Node::build(*node, Node::invalid_id);
141 98969 : new_node->processor_id() = elem->processor_id();
142 98969 : Node * added_node = mesh->add_node(std::move(new_node));
143 249021 : for (const auto j : elem->node_index_range())
144 249021 : if (elem->node_id(j) == node->id())
145 : {
146 98969 : elem->set_node(j, added_node);
147 98969 : break;
148 : }
149 :
150 : // Add boundary info to the new node
151 98969 : BoundaryInfo & boundary_info = mesh->get_boundary_info();
152 98969 : std::vector<boundary_id_type> node_boundary_ids;
153 98969 : boundary_info.boundary_ids(node, node_boundary_ids);
154 98969 : boundary_info.add_node(added_node, node_boundary_ids);
155 98969 : }
156 :
157 : void
158 49 : BreakMeshByElementGenerator::createInterface(MeshBase & mesh,
159 : const NodeToElemMapType & node_to_elem_map) const
160 : {
161 49 : std::set<std::pair<dof_id_type, unsigned int>> sides_breaking;
162 :
163 8740 : for (const auto & node_to_elems : node_to_elem_map)
164 116351 : for (const auto & elem_id_i : node_to_elems.second)
165 : {
166 107660 : Elem * elem_i = mesh.elem_ptr(elem_id_i);
167 1986420 : for (const auto & elem_id_j : node_to_elems.second)
168 : {
169 1878760 : Elem * elem_j = mesh.elem_ptr(elem_id_j);
170 1878760 : if (elem_i != elem_j && elem_i->has_neighbor(elem_j))
171 290384 : sides_breaking.insert(std::make_pair(elem_id_i, elem_i->which_neighbor_am_i(elem_j)));
172 : }
173 : }
174 :
175 49 : if (_interface_sides)
176 : {
177 49 : BoundaryInfo & boundary_info = mesh.get_boundary_info();
178 :
179 49 : const auto & existing_boundary_ids = boundary_info.get_boundary_ids();
180 : const boundary_id_type interface_id =
181 49 : existing_boundary_ids.empty() ? 0 : *existing_boundary_ids.rbegin() + 1;
182 49 : boundary_info.sideset_name(interface_id) = _interface_name;
183 :
184 97941 : for (const auto & [elem_id, side] : sides_breaking)
185 97892 : if (_interface_sides > 1 || elem_id > mesh.elem_ptr(elem_id)->neighbor_ptr(side)->id())
186 48946 : boundary_info.add_side(elem_id, side, interface_id);
187 : }
188 :
189 : // Remove element neighbor connections on the broken sides
190 97941 : for (const auto & [elem_id, side] : sides_breaking)
191 97892 : mesh.elem_ref(elem_id).set_neighbor(side, nullptr);
192 49 : }
|
