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