72{
73 auto const n_fractures = vec_fracture_mat_IDs.size();
74 std::map<unsigned, unsigned> matID_to_fid;
75 for (unsigned i = 0; i < n_fractures; i++)
76 {
77 matID_to_fid[vec_fracture_mat_IDs[i]] = i;
78 }
79
80
81 std::vector<std::size_t> all_fracture_nodes;
82 for (auto& vec : vec_fracture_nodes)
83 {
84 transform(cbegin(vec), cend(vec), back_inserter(all_fracture_nodes),
85 [](auto* const n) { return n->getID(); });
86 }
87
88
89 std::map<std::size_t, std::vector<std::size_t>> frac_nodeID_to_matIDs;
90 for (unsigned i = 0; i < n_fractures; i++)
91 {
92 for (auto* node : vec_fracture_nodes[i])
93 {
94 frac_nodeID_to_matIDs[node->getID()].push_back(
95 vec_fracture_mat_IDs[i]);
96 }
97 }
98
100
101
102 intersected_fracture_elements.resize(n_fractures);
103 for (auto frac_nodeID_to_matID : frac_nodeID_to_matIDs)
104 {
105 auto nodeID = frac_nodeID_to_matID.first;
106 auto const* node = mesh.getNode(frac_nodeID_to_matID.first);
107 auto const& matIDs = frac_nodeID_to_matID.second;
108 if (matIDs.size() < 2)
109 {
110 continue;
111 }
112
113 auto const elements_connected_to_node =
114 mesh.getElementsConnectedToNode(*node);
115 std::vector<MeshLib::Element*> conn_fracture_elements;
116 {
117 for (auto const* e : elements_connected_to_node)
118 {
119 if (e->getDimension() == (mesh.getDimension() - 1))
120 {
121 conn_fracture_elements.push_back(
122 const_cast<MeshLib::Element*>(e));
123 }
124 }
125 }
126
127 std::map<int, int> vec_matID_counts;
128 {
129 for (auto matid : matIDs)
130 {
131 vec_matID_counts[matid] = 0;
132 }
133
134 for (auto const* e : conn_fracture_elements)
135 {
136 auto matid = (*opt_material_ids)[e->getID()];
137 vec_matID_counts[matid]++;
138 }
139 }
140
141 for (auto matid : matIDs)
142 {
143 auto fid = matID_to_fid[matid];
144 for (auto* e : conn_fracture_elements)
145 {
146 auto e_matid = (*opt_material_ids)[e->getID()];
147 if (matID_to_fid[e_matid] != fid)
148 {
149 intersected_fracture_elements[fid].push_back(e);
150 }
151 }
152 }
153
154 bool isBranch = false;
155 {
156 for (auto vec_matID_count : vec_matID_counts)
157 {
158 auto count = vec_matID_count.second;
159 if (count % 2 == 1)
160 {
161 isBranch = true;
162 break;
163 }
164 }
165 }
166
167 if (isBranch)
168 {
169 std::vector<int> branch_matIDs(2);
170 for (auto vec_matID_count : vec_matID_counts)
171 {
172 auto matid = vec_matID_count.first;
173 auto count = vec_matID_count.second;
174 if (count % 2 == 0)
175 {
176 branch_matIDs[0] = matid;
177 }
178 else
179 {
180 branch_matIDs[1] = matid;
181 }
182 }
183 vec_branch_nodeID_matIDs.emplace_back(nodeID, branch_matIDs);
184 }
185 else
186 {
187 std::vector<int> junction_matIDs(2);
188 junction_matIDs[0] = std::min(vec_matID_counts.begin()->first,
189 vec_matID_counts.rbegin()->first);
190 junction_matIDs[1] = std::max(vec_matID_counts.begin()->first,
191 vec_matID_counts.rbegin()->first);
192 vec_junction_nodeID_matIDs.emplace_back(nodeID, junction_matIDs);
193 }
194 }
195
196 for (auto& elements : intersected_fracture_elements)
197 {
199 }
200
201 DBUG(
"-> found {:d} branches and {:d} junctions",
202 vec_branch_nodeID_matIDs.size(), vec_junction_nodeID_matIDs.size());
203}
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
void makeVectorUnique(std::vector< T > &v)
PropertyVector< int > const * materialIDs(Mesh const &mesh)