245 {
246
247 ParallelComm *pcomm =
249
250 MOFEM_LOG(
"EP", Sev::noisy) <<
"get_two_sides_of_crack_surface";
251
252 if (!pcomm->rank()) {
253
254 auto impl = [&](auto &saids) {
256
258
259 auto get_adj = [&](auto e, auto dim) {
262 e, dim, true, adj, moab::Interface::UNION),
263 "get adj");
264 return adj;
265 };
266
267 auto get_conn = [&](auto e) {
270 "get connectivity");
271 return conn;
272 };
273
274 constexpr bool debug =
false;
277 body_ents);
278 auto body_skin =
get_skin(m_field, body_ents);
279 auto body_skin_edges = get_adj(body_skin, 1);
280
281 auto crack_skin =
282 subtract(
get_skin(m_field, crack_faces), body_skin_edges);
283 auto crack_skin_conn = get_conn(crack_skin);
284 auto crack_skin_conn_edges = get_adj(crack_skin_conn, 1);
285 auto crack_edges = get_adj(crack_faces, 1);
286 crack_edges = subtract(crack_edges, crack_skin);
287 auto all_tets = get_adj(crack_edges, 3);
288 crack_edges = subtract(crack_edges, crack_skin_conn_edges);
289 auto crack_conn = get_conn(crack_edges);
290 all_tets.merge(get_adj(crack_conn, 3));
291
296 crack_edges);
297 }
298
299 if (crack_faces.size()) {
300 auto grow = [&](
auto r) {
301 auto crack_faces_conn = get_conn(crack_faces);
303 auto size_r = 0;
304 while (size_r !=
r.size() &&
r.size() > 0) {
306 CHKERR moab.get_connectivity(r,
v,
true);
307 v = subtract(
v, crack_faces_conn);
310 moab::Interface::UNION);
311 r = intersect(r, all_tets);
312 }
314 break;
315 }
316 }
318 };
319
320 Range all_tets_ord = all_tets;
321 while (all_tets.size()) {
322 Range faces = get_adj(unite(saids.first, saids.second), 2);
323 faces = subtract(crack_faces, faces);
324 if (faces.size()) {
326 auto fit = faces.begin();
327 for (; fit != faces.end(); ++fit) {
328 tets = intersect(get_adj(
Range(*fit, *fit), 3), all_tets);
329 if (tets.size() == 2) {
330 break;
331 }
332 }
333 if (tets.empty()) {
334 break;
335 } else {
336 saids.first.insert(tets[0]);
337 saids.first = grow(saids.first);
338 all_tets = subtract(all_tets, saids.first);
339 if (tets.size() == 2) {
340 saids.second.insert(tets[1]);
341 saids.second = grow(saids.second);
342 all_tets = subtract(all_tets, saids.second);
343 }
344 }
345 } else {
346 break;
347 }
348 }
349
350 saids.first = subtract(all_tets_ord, saids.second);
351 saids.second = subtract(all_tets_ord, saids.first);
352 }
353
355 };
356
357 std::pair<Range, Range> saids;
358 if (crack_faces.size())
360 return saids;
361 }
362
363 MOFEM_LOG(
"EP", Sev::noisy) <<
"get_two_sides_of_crack_surface <- done";
364
365 return std::pair<Range, Range>();
366}
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MOFEM_LOG(channel, severity)
Log.
const double v
phase velocity of light in medium (cm/ns)