259 {
260
261 ParallelComm *pcomm =
263
264 MOFEM_LOG(
"EP", Sev::noisy) <<
"get_two_sides_of_crack_surface";
265
266 if (!pcomm->rank()) {
267
268 auto impl = [&](auto &saids) {
270
272
273 auto get_adj = [&](auto e, auto dim) {
276 e, dim, true, adj, moab::Interface::UNION),
277 "get adj");
278 return adj;
279 };
280
281 auto get_conn = [&](auto e) {
284 "get connectivity");
285 return conn;
286 };
287
288 constexpr bool debug =
false;
291 body_ents);
292 auto body_skin =
get_skin(m_field, body_ents);
293 auto body_skin_edges = get_adj(body_skin, 1);
294
295 auto crack_skin =
296 subtract(
get_skin(m_field, crack_faces), body_skin_edges);
297 auto crack_skin_conn = get_conn(crack_skin);
298 auto crack_skin_conn_edges = get_adj(crack_skin_conn, 1);
299 auto crack_edges = get_adj(crack_faces, 1);
300 crack_edges = subtract(crack_edges, crack_skin);
301 auto all_tets = get_adj(crack_edges, 3);
302 crack_edges = subtract(crack_edges, crack_skin_conn_edges);
303 auto crack_conn = get_conn(crack_edges);
304 all_tets.merge(get_adj(crack_conn, 3));
305
310 crack_edges);
311 }
312
313 if (crack_faces.size()) {
314 auto grow = [&](
auto r) {
315 auto crack_faces_conn = get_conn(crack_faces);
317 auto size_r = 0;
318 while (size_r !=
r.size() &&
r.size() > 0) {
320 CHKERR moab.get_connectivity(r,
v,
true);
321 v = subtract(
v, crack_faces_conn);
324 moab::Interface::UNION);
325 r = intersect(r, all_tets);
326 }
328 break;
329 }
330 }
332 };
333
334 Range all_tets_ord = all_tets;
335 while (all_tets.size()) {
336 Range faces = get_adj(unite(saids.first, saids.second), 2);
337 faces = subtract(crack_faces, faces);
338 if (faces.size()) {
340 auto fit = faces.begin();
341 for (; fit != faces.end(); ++fit) {
342 tets = intersect(get_adj(
Range(*fit, *fit), 3), all_tets);
343 if (tets.size() == 2) {
344 break;
345 }
346 }
347 if (tets.empty()) {
348 break;
349 } else {
350 saids.first.insert(tets[0]);
351 saids.first = grow(saids.first);
352 all_tets = subtract(all_tets, saids.first);
353 if (tets.size() == 2) {
354 saids.second.insert(tets[1]);
355 saids.second = grow(saids.second);
356 all_tets = subtract(all_tets, saids.second);
357 }
358 }
359 } else {
360 break;
361 }
362 }
363
364 saids.first = subtract(all_tets_ord, saids.second);
365 saids.second = subtract(all_tets_ord, saids.first);
366 }
367
369 };
370
371 std::pair<Range, Range> saids;
372 if (crack_faces.size())
374 return saids;
375 }
376
377 MOFEM_LOG(
"EP", Sev::noisy) <<
"get_two_sides_of_crack_surface <- done";
378
379 return std::pair<Range, Range>();
380}
#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)