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 29856 : BreakMeshByElementGenerator::validParams()
24 : {
25 29856 : InputParameters params = MeshGenerator::validParams();
26 59712 : params.addClassDescription("Break all element-element interfaces in the specified subdomains.");
27 119424 : params.addRequiredParam<MeshGeneratorName>("input", "The mesh we want to modify");
28 89568 : params.addParam<std::vector<SubdomainID>>(
29 : "subdomains",
30 59712 : std::vector<SubdomainID>(),
31 : "The list of subdomain IDs to explode. Leave unset to explode all subdomains.");
32 119424 : params.addParam<BoundaryName>(
33 : "interface_name",
34 : "element_boundaries",
35 : "The boundary name containing all broken element-element interfaces.");
36 119424 : params.addRangeCheckedParam<unsigned int>(
37 : "interface_sides",
38 59712 : 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 29856 : return params;
43 0 : }
44 :
45 61 : BreakMeshByElementGenerator::BreakMeshByElementGenerator(const InputParameters & parameters)
46 : : MeshGenerator(parameters),
47 61 : _input(getMesh("input")),
48 122 : _subdomains(getParam<std::vector<SubdomainID>>("subdomains")),
49 122 : _interface_name(getParam<BoundaryName>("interface_name")),
50 183 : _interface_sides(getParam<unsigned int>("interface_sides"))
51 : {
52 61 : }
53 :
54 : std::unique_ptr<MeshBase>
55 58 : BreakMeshByElementGenerator::generate()
56 : {
57 58 : std::unique_ptr<MeshBase> mesh = std::move(_input);
58 :
59 58 : if (!mesh->is_prepared())
60 22 : mesh->prepare_for_use();
61 :
62 : // check that the subdomain IDs exist in the mesh
63 118 : for (const auto & id : _subdomains)
64 64 : if (!MooseMeshUtils::hasSubdomainID(*mesh, id))
65 8 : paramError("subdomains", "The block ID '", id, "' was not found in the mesh");
66 :
67 54 : BoundaryInfo & boundary_info = mesh->get_boundary_info();
68 108 : if (_interface_sides &&
69 54 : 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 54 : const auto node_to_elem_map = buildSubdomainRestrictedNodeToElemMap(mesh, _subdomains);
73 :
74 54 : duplicateNodes(mesh, node_to_elem_map);
75 :
76 54 : createInterface(*mesh, node_to_elem_map);
77 :
78 54 : 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 54 : boundary_info.regenerate_id_sets();
83 :
84 108 : return dynamic_pointer_cast<MeshBase>(mesh);
85 54 : }
86 :
87 : BreakMeshByElementGenerator::NodeToElemMapType
88 54 : BreakMeshByElementGenerator::buildSubdomainRestrictedNodeToElemMap(
89 : std::unique_ptr<MeshBase> & mesh, const std::vector<SubdomainID> & subdomains) const
90 : {
91 54 : NodeToElemMapType node_to_elem_map;
92 59118 : for (const auto & elem : mesh->active_element_ptr_range())
93 : {
94 : // Skip if subdomains are specified and the element is not in them
95 31032 : if (!subdomains.empty() &&
96 31032 : std::find(subdomains.begin(), subdomains.end(), elem->subdomain_id()) == subdomains.end())
97 216 : continue;
98 :
99 29316 : std::set<const Elem *> neighbors;
100 29316 : elem->find_point_neighbors(neighbors);
101 :
102 147156 : 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 117840 : bool should_duplicate = true;
107 117840 : if (!subdomains.empty())
108 47518 : for (auto neighbor : neighbors)
109 61119 : if (neighbor->contains_point(elem->node_ref(n)) &&
110 18953 : std::find(subdomains.begin(), subdomains.end(), neighbor->subdomain_id()) ==
111 61119 : subdomains.end())
112 : {
113 360 : should_duplicate = false;
114 360 : break;
115 : }
116 :
117 117840 : if (should_duplicate)
118 117480 : node_to_elem_map[elem->node_id(n)].insert(elem->id());
119 : }
120 29370 : }
121 :
122 54 : return node_to_elem_map;
123 0 : }
124 :
125 : void
126 54 : BreakMeshByElementGenerator::duplicateNodes(std::unique_ptr<MeshBase> & mesh,
127 : const NodeToElemMapType & node_to_elem_map) const
128 : {
129 9546 : for (const auto & [node_id, connected_elem_ids] : node_to_elem_map)
130 126972 : for (auto & connected_elem_id : connected_elem_ids)
131 117480 : if (connected_elem_id != *connected_elem_ids.begin())
132 107988 : duplicateNode(mesh, mesh->elem_ptr(connected_elem_id), mesh->node_ptr(node_id));
133 54 : }
134 :
135 : void
136 107988 : BreakMeshByElementGenerator::duplicateNode(std::unique_ptr<MeshBase> & mesh,
137 : Elem * elem,
138 : const Node * node) const
139 : {
140 107988 : std::unique_ptr<Node> new_node = Node::build(*node, Node::invalid_id);
141 107988 : new_node->processor_id() = elem->processor_id();
142 107988 : Node * added_node = mesh->add_node(std::move(new_node));
143 271740 : for (const auto j : elem->node_index_range())
144 271740 : if (elem->node_id(j) == node->id())
145 : {
146 107988 : elem->set_node(j, added_node);
147 107988 : break;
148 : }
149 :
150 : // Add boundary info to the new node
151 107988 : BoundaryInfo & boundary_info = mesh->get_boundary_info();
152 107988 : std::vector<boundary_id_type> node_boundary_ids;
153 107988 : boundary_info.boundary_ids(node, node_boundary_ids);
154 107988 : boundary_info.add_node(added_node, node_boundary_ids);
155 107988 : }
156 :
157 : void
158 54 : BreakMeshByElementGenerator::createInterface(MeshBase & mesh,
159 : const NodeToElemMapType & node_to_elem_map) const
160 : {
161 54 : std::set<std::pair<dof_id_type, unsigned int>> sides_breaking;
162 :
163 9546 : for (const auto & node_to_elems : node_to_elem_map)
164 126972 : for (const auto & elem_id_i : node_to_elems.second)
165 : {
166 117480 : Elem * elem_i = mesh.elem_ptr(elem_id_i);
167 2167170 : for (const auto & elem_id_j : node_to_elems.second)
168 : {
169 2049690 : Elem * elem_j = mesh.elem_ptr(elem_id_j);
170 2049690 : if (elem_i != elem_j && elem_i->has_neighbor(elem_j))
171 316842 : sides_breaking.insert(std::make_pair(elem_id_i, elem_i->which_neighbor_am_i(elem_j)));
172 : }
173 : }
174 :
175 54 : if (_interface_sides)
176 : {
177 54 : BoundaryInfo & boundary_info = mesh.get_boundary_info();
178 :
179 54 : const auto & existing_boundary_ids = boundary_info.get_boundary_ids();
180 : const boundary_id_type interface_id =
181 54 : existing_boundary_ids.empty() ? 0 : *existing_boundary_ids.rbegin() + 1;
182 54 : boundary_info.sideset_name(interface_id) = _interface_name;
183 :
184 106866 : for (const auto & [elem_id, side] : sides_breaking)
185 106812 : if (_interface_sides > 1 || elem_id > mesh.elem_ptr(elem_id)->neighbor_ptr(side)->id())
186 53406 : boundary_info.add_side(elem_id, side, interface_id);
187 : }
188 :
189 : // Remove element neighbor connections on the broken sides
190 106866 : for (const auto & [elem_id, side] : sides_breaking)
191 106812 : mesh.elem_ref(elem_id).set_neighbor(side, nullptr);
192 54 : }
|
