242 {
243
244 ParallelComm *pcomm =
246
247 MOFEM_LOG(
"EP", Sev::noisy) <<
"get_two_sides_of_crack_surface";
248
249 if (!pcomm->rank()) {
250
251 auto impl = [&](auto &saids) {
253
255
256 auto get_adj = [&](auto e, auto dim) {
259 e, dim, true, adj, moab::Interface::UNION),
260 "get adj");
261 return adj;
262 };
263
264 auto get_conn = [&](auto e) {
267 "get connectivity");
268 return conn;
269 };
270
271 constexpr bool debug =
false;
274 body_ents);
275 auto body_skin =
get_skin(m_field, body_ents);
276 auto body_skin_edges = get_adj(body_skin, 1);
277
278 auto crack_skin =
279 subtract(
get_skin(m_field, crack_faces), body_skin_edges);
280 auto crack_skin_conn = get_conn(crack_skin);
281 auto crack_skin_conn_edges = get_adj(crack_skin_conn, 1);
282 auto crack_edges = get_adj(crack_faces, 1);
283 crack_edges = subtract(crack_edges, crack_skin);
284 auto all_tets = get_adj(crack_edges, 3);
285 crack_edges = subtract(crack_edges, crack_skin_conn_edges);
286 auto crack_conn = get_conn(crack_edges);
287 all_tets.merge(get_adj(crack_conn, 3));
288
293 crack_edges);
294 }
295
296 if (crack_faces.size()) {
297 auto grow = [&](
auto r) {
298 auto crack_faces_conn = get_conn(crack_faces);
300 auto size_r = 0;
301 while (size_r !=
r.size() &&
r.size() > 0) {
303 CHKERR moab.get_connectivity(r,
v,
true);
304 v = subtract(
v, crack_faces_conn);
307 moab::Interface::UNION);
308 r = intersect(r, all_tets);
309 }
311 break;
312 }
313 }
315 };
316
317 Range all_tets_ord = all_tets;
318 while (all_tets.size()) {
319 Range faces = get_adj(unite(saids.first, saids.second), 2);
320 faces = subtract(crack_faces, faces);
321 if (faces.size()) {
323 auto fit = faces.begin();
324 for (; fit != faces.end(); ++fit) {
325 tets = intersect(get_adj(
Range(*fit, *fit), 3), all_tets);
326 if (tets.size() == 2) {
327 break;
328 }
329 }
330 if (tets.empty()) {
331 break;
332 } else {
333 saids.first.insert(tets[0]);
334 saids.first = grow(saids.first);
335 all_tets = subtract(all_tets, saids.first);
336 if (tets.size() == 2) {
337 saids.second.insert(tets[1]);
338 saids.second = grow(saids.second);
339 all_tets = subtract(all_tets, saids.second);
340 }
341 }
342 } else {
343 break;
344 }
345 }
346
347 saids.first = subtract(all_tets_ord, saids.second);
348 saids.second = subtract(all_tets_ord, saids.first);
349 }
350
352 };
353
354 std::pair<Range, Range> saids;
355 if (crack_faces.size())
357 return saids;
358 }
359
360 MOFEM_LOG(
"EP", Sev::noisy) <<
"get_two_sides_of_crack_surface <- done";
361
362 return std::pair<Range, Range>();
363}
#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)