24 #include "libmesh/int_range.h" 29 _faces(_num_faces, nullptr),
30 _face_neighbors(_num_faces),
31 _face_edge_neighbors(_num_faces)
43 EFAError(
"In EFAelement3D the supported TET element types are TET4, TET10 and TET14");
55 EFAError(
"In EFAelement3D the supported HEX element types are HEX8, HEX20 and HEX27");
58 EFAError(
"In EFAelement3D the supported element types are TET4, TET10, TET14, HEX8, HEX20 and " 64 :
EFAElement(from_elem->_id, from_elem->_num_nodes),
65 _num_faces(from_elem->_num_faces),
66 _faces(_num_faces, nullptr),
67 _face_neighbors(_num_faces),
68 _face_edge_neighbors(_num_faces)
82 EFAError(
"In EFAelement3D ",
84 " the copy constructor must have from_elem w/ global nodes. node: ",
87 from_elem->
_nodes[i]->category());
93 for (
unsigned int i = 0; i < from_elem->
_fragments.size(); ++i)
107 EFAError(
"In EFAelement3D copy constructor this elem's nodes must be local");
118 EFAError(
"this EFAelement3D constructor only converts global nodes to local nodes");
123 for (
unsigned int i = 0; i <
_fragments.size(); ++i)
127 for (
unsigned int i = 0; i <
_faces.size(); ++i)
268 EFAError(
"EFAElement3D: number of faces should be either 4(TET) or 6(HEX).");
284 bool node_in_frag =
false;
306 for (
unsigned int i = 0; i <
_nodes.size(); ++i)
307 non_physical_nodes.insert(
_nodes[i]);
310 std::set<EFANode *>::iterator sit;
311 for (sit = non_physical_nodes.begin(); sit != non_physical_nodes.end();)
314 for (
unsigned int i = 0; i <
_fragments.size(); ++i)
318 non_physical_nodes.erase(sit++);
334 if (
_nodes[i] == old_node)
337 for (
unsigned int i = 0; i <
_fragments.size(); ++i)
340 for (
unsigned int i = 0; i <
_faces.size(); ++i)
343 if (
_parent && descend_to_parent)
362 for (
unsigned int i = 0; i <
_fragments.size(); ++i)
374 std::vector<EFANode *> & master_nodes,
375 std::vector<double> & master_weights)
const 378 master_nodes.clear();
379 master_weights.clear();
380 bool masters_found =
false;
385 masters_found =
_faces[i]->getMasterInfo(node, master_nodes, master_weights);
389 EFAError(
"In getMasterInfo: cannot find master nodes in element faces");
399 std::vector<double> xi_3d(3, -100.0);
400 for (
unsigned int j = 0;
j < 3; ++
j)
404 master_nodes.push_back(
_nodes[
j]);
411 EFAError(
"unknown 3D element");
412 master_weights.push_back(
weight);
414 masters_found =
true;
421 EFAError(
"In EFAelement3D::getMasterInfo, cannot find the given EFANode");
433 bool overlays =
false;
436 EFAError(
"failed to dynamic cast to other3d");
440 if (common_face_curr.size() == 1)
442 unsigned int curr_face_id = common_face_curr[0];
444 unsigned int other_face_id = other3d->
getFaceID(curr_face);
449 else if (common_face_curr.size() > 1)
464 EFAError(
"in getNeighborIndex() element ",
_id,
" does not have neighbor ", neighbor_elem->
id());
469 unsigned int face_id,
470 unsigned int edge_id,
471 unsigned int & neigh_face_id,
472 unsigned int & neigh_edge_id)
const 475 for (
unsigned int i = 0; i < neighbor_elem->
numFaces(); ++i)
488 EFAError(
"in getNeighborEdgeIndex() element ",
490 " does not share a common edge with element",
491 neighbor_elem->
id());
498 for (
unsigned int face_iter = 0; face_iter <
_num_faces; ++face_iter)
501 for (
unsigned int edge_iter = 0; edge_iter <
_faces[face_iter]->numEdges(); ++edge_iter)
510 std::vector<std::pair<unsigned int, unsigned int>> common_ids;
515 EFAError(
"neighbor_elem is not of EFAelement3D type");
518 if (common_face_id.empty() &&
getCommonEdgeID(neigh_elem, common_ids) &&
521 bool is_edge_neighbor =
false;
525 EFAError(
"in updateFaceNeighbors: Cannot have more than 1 fragment");
530 is_edge_neighbor =
true;
533 is_edge_neighbor =
true;
535 if (is_edge_neighbor)
536 for (
const auto & [face_id, edge_id] : common_ids)
542 unsigned int face_id = common_face_id[0];
543 bool is_face_neighbor =
false;
547 EFAError(
"in updateFaceNeighbors: Cannot have more than 1 fragment");
552 is_face_neighbor =
true;
555 is_face_neighbor =
true;
557 if (is_face_neighbor)
562 " already has 2 face neighbors: ",
575 for (
unsigned int face_iter = 0; face_iter <
_num_faces; ++face_iter)
577 for (
unsigned int en_iter = 0; en_iter <
_face_neighbors[face_iter].size(); ++en_iter)
580 if (neigh_elem !=
nullptr)
582 bool found_neighbor =
false;
583 for (
unsigned int face_iter2 = 0; face_iter2 < neigh_elem->
numFaces(); ++face_iter2)
585 for (
unsigned int en_iter2 = 0; en_iter2 < neigh_elem->
numFaceNeighbors(face_iter2);
590 if ((en_iter2 > 1) && (en_iter > 1))
592 "Element and neighbor element cannot both have >1 neighbors on a common face");
593 found_neighbor =
true;
599 EFAError(
"Neighbor element doesn't recognize current element as neighbor");
610 CrackTipElements.insert(
this);
611 for (
unsigned int face_iter = 0; face_iter <
_num_faces; ++face_iter)
620 EFAError(
"Element has a neighbor that overlays itself");
624 EFAError(
"crack_tip_split_element already flagged. In elem: ",
652 bool should_duplicate =
false;
655 std::set<EFAElement *>::iterator sit;
656 sit = CrackTipElements.find(
this);
658 should_duplicate =
true;
660 should_duplicate =
true;
662 should_duplicate =
true;
664 return should_duplicate;
674 bool should_duplicate =
false;
677 std::vector<unsigned int> split_neighbors;
679 should_duplicate =
true;
685 std::set<EFANode *> non_physical_nodes;
692 EFAError(
"general elem is not of type EFAelement3D");
695 std::set<EFAElement *>::iterator sit;
696 sit = CrackTipElements.find(neigh_elem);
697 if (sit != CrackTipElements.end() && neigh_elem->
numFragments() > 1)
699 for (
unsigned int i = 0; i < neigh_elem->
numFaces(); ++i)
701 std::set<EFANode *> neigh_face_nodes = neigh_elem->
getFaceNodes(i);
705 should_duplicate =
true;
710 if (should_duplicate)
715 return should_duplicate;
724 bool should_duplicate =
false;
735 for (
unsigned int j = 0;
j < neighbor_elem->
numFaces(); ++
j)
743 should_duplicate =
true;
750 if (should_duplicate)
754 return should_duplicate;
763 bool will_extend =
false;
770 EFAError(
"in will_crack_tip_extend() element ",
774 " neighbors on face ",
780 EFAError(
"in will_crack_tip_extend() element ",
790 std::vector<EFANode *> neigh_cut_nodes = neigh_frag1->
getCommonNodes(neigh_frag2);
791 unsigned int counter = 0;
792 for (
unsigned int j = 0;
j < neigh_cut_nodes.size(); ++
j)
799 split_neighbors.push_back(neigh_idx);
817 unsigned int num_cut_faces = 0;
819 if (
_faces[i]->hasIntersection())
821 return num_cut_faces;
829 bool cut_third =
false;
832 unsigned int num_interior_faces = 0;
833 for (
unsigned int i = 0; i <
_fragments[0]->numFaces(); ++i)
836 num_interior_faces += 1;
838 if (num_interior_faces == 3)
846 std::map<unsigned int, EFANode *> & EmbeddedNodes)
849 std::set<EFAElement *>::iterator sit;
850 sit = CrackTipElements.find(
this);
851 if (sit != CrackTipElements.end())
856 EFAError(
"crack tip elem ",
_id,
" must have 1 fragment");
862 _fragments[0]->removeInvalidEmbeddedNodes(EmbeddedNodes);
868 EFAError(
"Element ",
_id,
" must have 1 fragment at this point");
871 unsigned int num_cut_frag_faces =
_fragments[0]->getNumCuts();
872 unsigned int num_frag_faces =
_fragments[0]->numFaces();
874 EFAError(
"In element ",
_id,
" there are too many cut fragment faces");
877 if (num_cut_frag_faces == 0)
888 std::vector<EFAFragment3D *> new_frags =
_fragments[0]->split();
889 if (new_frags.size() == 1 || new_frags.size() == 2)
893 for (
unsigned int i = 0; i < new_frags.size(); ++i)
897 EFAError(
"Number of fragments must be 1 or 2 at this point");
906 if (n_interior_nodes > 0 && n_interior_nodes != 1)
907 EFAError(
"After update_fragments this element has ", n_interior_nodes,
" interior nodes");
909 if (n_old_frag_cuts == 0)
912 EFAError(
"Incorrect link size for element with 0 cuts");
917 EFAError(
"Incorrect link size for element with crack-tip faces");
922 n_old_frag_faces + n_old_frag_cuts + 2)
923 EFAError(
"Incorrect link size for element that has been completely cut");
932 EFAError(
"from_elem is not of EFAelement3D type");
936 EFAError(
"in restoreFragmentInfo elements must not have any pre-existing fragments");
937 for (
unsigned int i = 0; i < from_elem3d->
numFragments(); ++i)
942 EFAError(
"in restoreFragmentInfo elements must not have any pre-exsiting interior nodes");
948 EFAError(
"In restoreEdgeIntersection: edge cuts already exist in element ",
_id);
953 for (
unsigned int i = 0; i < from_elem3d->
numNodes(); ++i)
959 EFAError(
"In restoreFragmentInfo all of from_elem's nodes must be local");
965 std::map<unsigned int, EFAElement *> & Elements,
966 std::map<unsigned int, EFAElement *> & newChildElements,
967 std::vector<EFAElement *> & ChildElements,
968 std::vector<EFAElement *> & ParentElements,
969 std::map<unsigned int, EFANode *> & TempNodes)
972 EFAError(
"Element cannot have existing children in createChildElements");
977 EFAError(
"More than 2 fragments not yet supported");
980 ParentElements.push_back(
this);
981 for (
unsigned int ichild = 0; ichild <
_fragments.size(); ++ichild)
983 unsigned int new_elem_id;
984 if (newChildElements.size() == 0)
990 newChildElements.insert(std::make_pair(new_elem_id, childElem));
992 ChildElements.push_back(childElem);
996 std::vector<std::vector<EFANode *>> cut_plane_nodes;
1002 std::vector<EFANode *> node_line;
1003 for (
unsigned int j = 0;
j < face->
numNodes(); ++
j)
1004 node_line.push_back(face->
getNode(
j));
1005 cut_plane_nodes.push_back(node_line);
1009 std::vector<EFAPoint> cut_plane_points;
1014 if (cut_plane_nodes.size())
1016 for (
unsigned int i = 0; i < cut_plane_nodes[0].size(); ++i)
1018 std::vector<EFANode *> master_nodes;
1019 std::vector<double> master_weights;
1021 this->
getMasterInfo(cut_plane_nodes[0][i], master_nodes, master_weights);
1023 for (
unsigned int i = 0; i < master_nodes.size(); ++i)
1029 cut_plane_points.push_back(coor);
1031 for (
unsigned int i = 0; i < cut_plane_points.size(); ++i)
1032 orig += cut_plane_points[i];
1033 orig /= cut_plane_points.size();
1036 for (
unsigned int i = 0; i < cut_plane_points.size(); ++i)
1037 center += cut_plane_points[i];
1038 center /= cut_plane_points.size();
1040 for (
unsigned int i = 0; i < cut_plane_points.size(); ++i)
1042 unsigned int iplus1 = i < cut_plane_points.size() - 1 ? i + 1 : 0;
1045 normal += ray1.
cross(ray2);
1047 normal /= cut_plane_points.size();
1056 EFAPoint origin_to_point = p - orig;
1059 else if (origin_to_point * normal <
Xfem::tol)
1065 TempNodes.insert(std::make_pair(new_node_id, newNode));
1098 for (
unsigned int i = 0; i <
_num_faces; ++i)
1101 std::vector<EFANode *> nodes_to_delete;
1102 for (
unsigned int j = 0;
j <
_faces[i]->numEdges(); ++
j)
1105 for (
unsigned int k = 0;
k <
edge->numEmbeddedNodes(); ++
k)
1108 nodes_to_delete.push_back(
edge->getEmbeddedNode(
k));
1113 for (
unsigned int j = 0;
j <
_faces[i]->numInteriorNodes(); ++
j)
1117 nodes_to_delete.push_back(face_node);
1121 for (
unsigned int j = 0;
j < nodes_to_delete.size(); ++
j)
1129 std::map<unsigned int, EFANode *> & TempNodes,
1130 std::map<
EFANode *, std::set<EFAElement *>> & InverseConnectivityMap,
1131 bool merge_phantom_faces)
1135 EFAError(
"no parent element for child element ",
_id,
" in connect_neighbors");
1138 EFAError(
"cannot dynamic cast to parent3d in connect_neighbors");
1148 if (
_faces[
j]->hasIntersection())
1150 for (
unsigned int l = 0; l < NeighborElem->
numChildren(); ++l)
1154 if (!childOfNeighborElem)
1155 EFAError(
"dynamic cast childOfNeighborElem fails");
1158 EFAFace * neighborChildFace = childOfNeighborElem->
getFace(neighbor_face_id);
1159 if (
_faces[
j]->equivalent(neighborChildFace))
1164 for (
unsigned int i = 0; i <
_faces[
j]->numNodes(); ++i)
1166 unsigned int childNodeIndex = i;
1167 unsigned int neighborChildNodeIndex =
1171 EFANode * childOfNeighborNode = neighborChildFace->
getNode(neighborChildNodeIndex);
1173 childNode, childOfNeighborNode, childOfNeighborElem, PermanentNodes, TempNodes);
1178 unsigned int childNodeIndex = m;
1179 unsigned int neighborChildNodeIndex =
1182 EFANode * childNode =
_faces[
j]->getInteriorFaceNode(childNodeIndex);
1183 EFANode * childOfNeighborNode =
1186 childNode, childOfNeighborNode, childOfNeighborElem, PermanentNodes, TempNodes);
1194 if (merge_phantom_faces)
1196 for (
unsigned int l = 0; l < NeighborElem->
numChildren(); ++l)
1200 if (!childOfNeighborElem)
1201 EFAError(
"dynamic cast childOfNeighborElem fails");
1203 EFAFace * neighborChildFace = childOfNeighborElem->
getFace(neighbor_face_id);
1204 if (!neighborChildFace
1205 ->hasIntersection())
1209 if (
_faces[
j]->equivalent(neighborChildFace))
1212 for (
unsigned int i = 0; i <
_faces[
j]->numNodes(); ++i)
1214 unsigned int childNodeIndex = i;
1215 unsigned int neighborChildNodeIndex =
1219 EFANode * childOfNeighborNode = neighborChildFace->
getNode(neighborChildNodeIndex);
1221 if (childNode->
parent() !=
nullptr &&
1226 childOfNeighborNode,
1227 childOfNeighborElem,
1248 std::set<EFAElement *> patch_elems = InverseConnectivityMap[childNode->
parent()];
1249 if (parent3d->
numFragments() == 1 && patch_elems.size() == 1)
1256 PermanentNodes.insert(std::make_pair(new_node_id, newNode));
1261 "Attempted to delete node: ", childNode->
id(),
" from TempNodes, but couldn't find it");
1270 ostream << std::setw(5);
1271 ostream <<
_id <<
"| ";
1274 for (
unsigned int k = 0;
k <
_faces[
j]->numNodes(); ++
k)
1275 ostream << std::setw(5) <<
_faces[
j]->getNode(
k)->idCatString();
1278 ostream << std::endl;
1281 ostream << std::setw(5);
1286 for (
unsigned int k = 0;
k <
_faces[
j]->numEdges(); ++
k)
1288 ostream << std::setw(4);
1289 if (
_faces[
j]->getEdge(
k)->hasIntersection())
1291 if (
_faces[
j]->getEdge(
k)->numEmbeddedNodes() > 1)
1294 for (
unsigned int l = 0; l <
_faces[
j]->getEdge(
k)->numEmbeddedNodes(); ++l)
1296 ostream <<
_faces[
j]->getEdge(
k)->getEmbeddedNode(l)->id() <<
" ";
1297 if (l ==
_faces[
j]->getEdge(
k)->numEmbeddedNodes() - 1)
1304 ostream <<
_faces[
j]->getEdge(
k)->getEmbeddedNode(0)->id() <<
" ";
1311 ostream << std::endl;
1314 ostream << std::setw(5);
1319 ostream << std::setw(4);
1340 ostream << std::endl;
1345 ostream << std::setw(4);
1346 ostream <<
"frag" <<
j <<
"| ";
1349 for (
unsigned int l = 0; l <
_fragments[
j]->getFace(
k)->numNodes(); ++l)
1350 ostream << std::setw(5) <<
_fragments[
j]->getFace(
k)->getNode(l)->idCatString();
1353 ostream << std::endl;
1355 ostream << std::endl;
1364 EFAError(
"frag_id out of bounds");
1370 std::set<EFANode *> face_nodes;
1371 for (
unsigned int i = 0; i <
_faces[face_id]->numNodes(); ++i)
1372 face_nodes.insert(
_faces[face_id]->getNode(i));
1380 unsigned int face_id = std::numeric_limits<unsigned int>::max();
1381 bool face_found =
false;
1382 for (
unsigned int i = 0; i <
_num_faces; ++i)
1393 std::vector<double> xi_2d(2, 0.0);
1394 if (
_faces[face_id]->getFaceNodeParametricCoords(node, xi_2d))
1397 EFAError(
"failed to get the 2D para coords on the face");
1408 EFAError(
"interior_node_id out of bounds");
1414 for (
unsigned int i = 0; i <
_fragments.size(); ++i)
1417 for (
unsigned int i = 0; i <
_faces.size(); ++i)
1420 if (remove_for_neighbor)
1422 for (
unsigned int i = 0; i <
numFaces(); ++i)
1444 int hex_local_node_indices[6][4] = {
1445 {0, 3, 2, 1}, {0, 1, 5, 4}, {1, 2, 6, 5}, {2, 3, 7, 6}, {3, 0, 4, 7}, {4, 5, 6, 7}};
1446 int tet_local_node_indices[4][3] = {{0, 2, 1}, {0, 1, 3}, {1, 2, 3}, {2, 0, 3}};
1448 int hex_interior_face_node_indices[6][5] = {{8, 9, 10, 11, 20},
1449 {8, 13, 16, 12, 21},
1450 {9, 14, 17, 13, 22},
1451 {10, 14, 18, 15, 23},
1452 {11, 15, 19, 12, 24},
1453 {16, 17, 18, 19, 25}};
1454 int tet_interior_face_node_indices[4][4] = {
1455 {4, 5, 6, 10}, {4, 7, 8, 11}, {5, 8, 9, 12}, {6, 7, 9, 13}};
1461 EFAError(
"num_faces of hexes must be 6");
1462 for (
unsigned int i = 0; i <
_num_faces; ++i)
1465 for (
unsigned int j = 0;
j < 4; ++
j)
1467 _faces[i]->createEdges();
1469 _faces[i]->setInteriorFaceNode(
k,
_nodes[hex_interior_face_node_indices[i][
k]]);
1475 EFAError(
"num_faces of tets must be 4");
1476 for (
unsigned int i = 0; i <
_num_faces; ++i)
1479 for (
unsigned int j = 0;
j < 3; ++
j)
1481 _faces[i]->createEdges();
1483 _faces[i]->setInteriorFaceNode(
k,
_nodes[tet_interior_face_node_indices[i][
k]]);
1487 EFAError(
"unknown 3D element type in createFaces()");
1489 for (
unsigned int face_iter = 0; face_iter <
_num_faces; ++face_iter)
1492 for (
unsigned int edge_iter = 0; edge_iter <
_faces[face_iter]->numEdges(); ++edge_iter)
1509 bool found_face_id =
false;
1510 unsigned int face_id;
1511 for (
unsigned int iface = 0; iface <
_num_faces; ++iface)
1513 if (
_faces[iface]->equivalent(face))
1516 found_face_id =
true;
1521 EFAError(
"input face not found in get_face_id()");
1525 std::vector<unsigned int>
1528 std::vector<unsigned int> face_id;
1529 for (
unsigned int i = 0; i <
_num_faces; ++i)
1533 face_id.push_back(i);
1542 std::vector<std::pair<unsigned int, unsigned int>> & common_ids)
const 1545 std::set<std::pair<EFANode *, EFANode *>> other_edges;
1548 const auto & face = *other_elem->
_faces[
k];
1551 const auto &
edge = *face.getEdge(l);
1552 other_edges.insert(
edge.getSortedNodes());
1562 const auto edge_nodes =
edge.getSortedNodes();
1565 if (
edge.isEmbeddedPermanent() || other_edges.count(edge_nodes) == 0)
1568 common_ids.emplace_back(i,
j);
1571 return common_ids.size() > 0;
1576 unsigned int node_id,
1580 bool found_id =
false;
1581 unsigned int neigh_face_node_id;
1583 if (common_face_id == face_id)
1587 for (
unsigned int i = 0; i < neigh_face->
numNodes(); ++i)
1591 neigh_face_node_id = i;
1598 EFAError(
"getNeighborFaceNodeID: neighbor_elem is not a neighbor on face_id");
1600 EFAError(
"getNeighborFaceNodeID: could not find neighbor face node id");
1601 return neigh_face_node_id;
1606 unsigned int node_id,
1610 bool found_id =
false;
1611 unsigned int neigh_face_node_id;
1613 if (common_face_id == face_id)
1622 neigh_face_node_id = i;
1629 EFAError(
"getNeighborFaceNodeID: neighbor_elem is not a neighbor on face_id");
1631 EFAError(
"getNeighborFaceNodeID: could not find neighbor face node id");
1632 return neigh_face_node_id;
1637 unsigned int edge_id,
1641 bool found_id =
false;
1642 unsigned int neigh_face_edge_id;
1644 if (common_face_id == face_id)
1648 for (
unsigned int i = 0; i < neigh_face->
numEdges(); ++i)
1650 if (
_faces[face_id]->getEdge(edge_id)->equivalent(*neigh_face->
getEdge(i)))
1652 neigh_face_edge_id = i;
1659 EFAError(
"getNeighborFaceEdgeID: neighbor_elem is not a neighbor on face_id");
1661 EFAError(
"getNeighborFaceEdgeID: could not find neighbor face edge id");
1662 return neigh_face_edge_id;
1669 for (
unsigned int i = 0; i <
_faces.size(); ++i)
1671 std::vector<EFAFace *> face_adjacents(
_faces[i]->numEdges(),
nullptr);
1672 for (
unsigned int j = 0;
j <
_faces.size(); ++
j)
1676 unsigned int adj_edge =
_faces[i]->adjacentCommonEdge(
_faces[
j]);
1677 face_adjacents[adj_edge] =
_faces[
j];
1696 EFAError(
"frag_id out of bounds in getFragmentFace");
1702 std::set<EFANode *> phantom_nodes;
1705 for (
unsigned int j = 0;
j <
_faces[face_id]->numNodes(); ++
j)
1707 bool node_in_frag =
false;
1712 node_in_frag =
true;
1720 return phantom_nodes;
1729 bool frag_face_found =
false;
1730 frag_face_id = std::numeric_limits<unsigned int>::max();
1733 for (
unsigned int j = 0;
j <
_fragments[0]->numFaces(); ++
j)
1738 frag_face_found =
true;
1743 return frag_face_found;
1748 unsigned int ElemFaceEdgeID,
1749 unsigned int & FragFaceID,
1750 unsigned int & FragFaceEdgeID)
const 1753 bool frag_edge_found =
false;
1754 FragFaceID = FragFaceEdgeID = std::numeric_limits<unsigned int>::max();
1757 EFAEdge * elem_edge =
_faces[ElemFaceID]->getEdge(ElemFaceEdgeID);
1759 for (
unsigned int i = 0; i < frag_face->
numEdges(); ++i)
1764 frag_edge_found =
true;
1769 return frag_edge_found;
1774 unsigned int ElemFaceEdgeID,
1775 double position)
const 1777 unsigned int FragFaceID, FragFaceEdgeID = std::numeric_limits<unsigned int>::max();
1778 bool is_in_real =
false;
1785 EFAEdge * elem_edge =
_faces[ElemFaceID]->getEdge(ElemFaceEdgeID);
1787 double xi[2] = {-1.0, -1.0};
1790 if ((position - xi[0]) * (position - xi[1]) <
1800 bool is_phantom =
false;
1803 bool contains_frag_face =
false;
1804 for (
unsigned int i = 0; i <
_fragments.size(); ++i)
1806 for (
unsigned int j = 0;
j <
_fragments[i]->numFaces(); ++
j)
1810 contains_frag_face =
true;
1814 if (contains_frag_face)
1817 if (!contains_frag_face)
1841 EFAError(
"edge neighbor does not exist");
1846 unsigned int edge_id,
1847 unsigned int neighbor_id)
const 1852 EFAError(
"edge neighbor does not exist");
1858 bool has_tip_faces =
false;
1861 for (
unsigned int i = 0; i <
_num_faces; ++i)
1863 unsigned int num_frag_faces = 0;
1864 if (
_faces[i]->hasIntersection())
1866 for (
unsigned int j = 0;
j <
_fragments[0]->numFaces(); ++
j)
1869 num_frag_faces += 1;
1871 if (num_frag_faces == 2)
1873 has_tip_faces =
true;
1879 return has_tip_faces;
1882 std::vector<unsigned int>
1886 std::vector<unsigned int> tip_face_id;
1889 for (
unsigned int i = 0; i <
_num_faces; ++i)
1891 unsigned int num_frag_faces = 0;
1892 if (
_faces[i]->hasIntersection())
1894 for (
unsigned int j = 0;
j <
_fragments[0]->numFaces(); ++
j)
1897 num_frag_faces += 1;
1899 if (num_frag_faces == 2)
1900 tip_face_id.push_back(i);
1911 std::set<EFANode *> tip_emb;
1914 for (
unsigned int i = 0; i <
_num_faces; ++i)
1916 std::vector<EFAFace *> frag_faces;
1917 if (
_faces[i]->hasIntersection())
1919 for (
unsigned int j = 0;
j <
_fragments[0]->numFaces(); ++
j)
1924 if (frag_faces.size() == 2)
1926 unsigned int edge_id = frag_faces[0]->adjacentCommonEdge(frag_faces[1]);
1927 tip_emb.insert(frag_faces[0]->getEdge(edge_id)->
getNode(0));
1928 tip_emb.insert(frag_faces[0]->getEdge(edge_id)->
getNode(1));
1939 bool contain_tip =
false;
1942 unsigned int num_frag_faces = 0;
1943 if (
_faces[face_id]->hasIntersection())
1945 for (
unsigned int j = 0;
j <
_fragments[0]->numFaces(); ++
j)
1948 num_frag_faces += 1;
1950 if (num_frag_faces == 2)
1961 bool has_cut =
false;
1966 unsigned int FragFaceID = std::numeric_limits<unsigned int>::max();
1979 unsigned int edge_id,
1982 std::map<unsigned int, EFANode *> & EmbeddedNodes,
1983 bool add_to_neighbor,
1984 bool add_to_adjacent)
1987 EFANode * local_embedded =
nullptr;
1991 local_embedded = embedded_node;
1995 unsigned int adj_face_id =
getFaceID(adj_face);
1999 bool cut_exist =
false;
2005 if (embedded_node && embedded_node != old_emb)
2006 EFAError(
"Attempting to add edge intersection when one already exists with different node.",
2013 local_embedded = old_emb;
2023 adj_face_id, adj_edge_id, position, edge_node1, embedded_node, local_embedded);
2025 if (!local_embedded)
2029 EmbeddedNodes.insert(std::make_pair(new_node_id, local_embedded));
2035 EFAError(
"element edge can't have >2 embedded nodes");
2039 unsigned int FragFaceID;
2040 unsigned int FragFaceEdgeID;
2045 double xi[2] = {-1.0, -1.0};
2048 double frag_pos = (position - xi[0]) / (xi[1] - xi[0]);
2052 frag_edge->
addIntersection(frag_pos, local_embedded, frag_edge_node1);
2056 if (add_to_adjacent)
2058 double adj_pos = 1.0 - position;
2060 adj_face_id, adj_edge_id, adj_pos, local_embedded, EmbeddedNodes,
false,
false);
2065 if (add_to_neighbor)
2067 for (
unsigned int en_iter = 0; en_iter <
numFaceNeighbors(face_id); ++en_iter)
2072 double neigh_pos = 1.0 - position;
2074 neigh_face_id, neigh_edge_id, neigh_pos, local_embedded, EmbeddedNodes,
false,
true);
2077 for (
unsigned int en_iter = 0; en_iter <
numEdgeNeighbors(face_id, edge_id); ++en_iter)
2080 unsigned int neigh_face_id, neigh_edge_id;
2087 neigh_pos = position;
2088 else if (
_faces[face_id]->getEdge(edge_id)->getNode(1) ==
2090 neigh_pos = 1.0 - position;
2092 EFAError(
"The EFANodes on commaon edge are not matched.");
2095 neigh_face_id, neigh_edge_id, neigh_pos, local_embedded, EmbeddedNodes,
false,
true);
2104 std::map<unsigned int, EFANode *> & ,
2115 unsigned int edge_id,
2125 for (
unsigned int en_iter = 0; en_iter <
numFaceNeighbors(face_id); ++en_iter)
2137 if (embedded_node && embedded_node != old_emb)
2139 "attempting to add edge intersection when one already exists with different node.");
2140 if (local_embedded && local_embedded != old_emb)
2141 EFAError(
"attempting to assign contradictory pointer to local_embedded.");
2143 local_embedded = old_emb;
2150 std::vector<double> & xi_2d,
2151 std::vector<double> & xi_3d)
const 2154 xi_3d.resize(3, 0.0);
2159 xi_3d[0] = xi_2d[1];
2160 xi_3d[1] = xi_2d[0];
2163 else if (face_id == 1)
2165 xi_3d[0] = xi_2d[0];
2167 xi_3d[2] = xi_2d[1];
2169 else if (face_id == 2)
2172 xi_3d[1] = xi_2d[0];
2173 xi_3d[2] = xi_2d[1];
2175 else if (face_id == 3)
2177 xi_3d[0] = -xi_2d[0];
2179 xi_3d[2] = xi_2d[1];
2181 else if (face_id == 4)
2184 xi_3d[1] = -xi_2d[0];
2185 xi_3d[2] = xi_2d[1];
2187 else if (face_id == 5)
2189 xi_3d[0] = xi_2d[0];
2190 xi_3d[1] = xi_2d[1];
2194 EFAError(
"face_id out of bounds");
2200 xi_3d[0] = xi_2d[0];
2201 xi_3d[1] = xi_2d[1];
2204 else if (face_id == 1)
2207 xi_3d[1] = xi_2d[0];
2208 xi_3d[2] = xi_2d[1];
2210 else if (face_id == 2)
2212 xi_3d[0] = xi_2d[1];
2214 xi_3d[2] = xi_2d[0];
2216 else if (face_id == 3)
2218 xi_3d[0] = xi_2d[0];
2219 xi_3d[1] = xi_2d[2];
2220 xi_3d[2] = xi_2d[1];
2223 EFAError(
"face_id out of bounds");
2226 EFAError(
"unknown element for 3D");
2229 std::vector<EFANode *>
2232 std::set<EFANode *> e1nodes(
_nodes.begin(),
2234 std::set<EFANode *> e2nodes(other_elem->
_nodes.begin(),
2237 return common_nodes;
void mapParametricCoordinateFrom2DTo3D(unsigned int face_id, std::vector< double > &xi_2d, std::vector< double > &xi_3d) const
EFAFragment3D * getFragment(unsigned int frag_id) const
virtual void connectNeighbors(std::map< unsigned int, EFANode *> &PermanentNodes, std::map< unsigned int, EFANode *> &TempNodes, std::map< EFANode *, std::set< EFAElement *>> &InverseConnectivityMap, bool merge_phantom_faces)
std::vector< EFAFace * > _faces
virtual bool shouldDuplicateForCrackTip(const std::set< EFAElement *> &CrackTipElements)
bool isCrackTipSplit() const
std::vector< std::vector< std::vector< EFAElement3D * > > > _face_edge_neighbors
EFAFace * getFace(unsigned int face_id) const
EFAElement3D * getFaceNeighbor(unsigned int face_id, unsigned int neighbor_id) const
EFANode * getEmbeddedNode(unsigned int index) const
virtual void switchNode(EFANode *new_node, EFANode *old_node, bool descend_to_parent)
std::vector< std::vector< EFAFace * > > _faces_adjacent_to_faces
unsigned int numChildren() const
unsigned int numEdges() const
bool isFacePhantom(unsigned int face_id) const
unsigned int getNewID(std::map< unsigned int, T *> &theMap)
bool hasIntersection() const
virtual unsigned int numInteriorNodes() const
virtual void printElement(std::ostream &ostream) const
bool getFragmentFaceEdgeID(unsigned int ElemFaceID, unsigned int ElemFaceEdgeID, unsigned int &FragFaceID, unsigned int &FragFaceEdgeID) const
void checkNeighborFaceCut(unsigned int face_id, unsigned int edge_id, double position, EFANode *from_node, EFANode *embedded_node, EFANode *&local_embedded)
bool getFaceNodeParametricCoordinates(EFANode *node, std::vector< double > &xi_3d) const
virtual bool shouldDuplicateForPhantomCorner()
bool getCommonEdgeID(const EFAElement3D *other_elem, std::vector< std::pair< unsigned int, unsigned int >> &common_ids) const
unsigned int _num_interior_face_nodes
void findFacesAdjacentToFaces()
std::vector< std::vector< EFAElement3D * > > _face_neighbors
bool hasSameOrientation(const EFAFace *other_face) const
unsigned int getFaceID(EFAFace *face) const
EFAFace * getAdjacentFace(unsigned int face_id, unsigned int edge_id) const
virtual void getNonPhysicalNodes(std::set< EFANode *> &non_physical_nodes) const
bool hasIntersection() const
bool deleteFromMap(std::map< unsigned int, T *> &theMap, T *elemToDelete, bool delete_elem=true)
virtual unsigned int numFragments() const
EFAPoint cross(const EFAPoint &point)
EFAElement3D * getEdgeNeighbor(unsigned int face_id, unsigned int edge_id, unsigned int neighbor_id) const
bool equivalent(const EFAFace *other_face) const
unsigned int getNeighborFaceNodeID(unsigned int face_id, unsigned int node_id, EFAElement3D *neighbor_elem) const
virtual void updateFragments(const std::set< EFAElement *> &CrackTipElements, std::map< unsigned int, EFANode *> &EmbeddedNodes)
unsigned int numGeneralNeighbors() const
unsigned int numFaces() const
void setFace(unsigned int face_id, EFAFace *face)
EFANode * getNode(unsigned int node_id) const
bool containsNode(EFANode *node) const
virtual unsigned int getNeighborIndex(const EFAElement *neighbor_elem) const
virtual void initCrackTip(std::set< EFAElement *> &CrackTipElements)
unsigned int getNeighborFaceEdgeID(unsigned int face_id, unsigned int edg_id, EFAElement3D *neighbor_elem) const
unsigned int numEdgeNeighbors(unsigned int face_id, unsigned int edge_id) const
EFAElement * getChild(unsigned int child_id) const
double distanceFromNode1(EFANode *node) const
std::vector< EFAPoint > _local_node_coor
std::vector< EFANode * > getCommonNodes(EFAFragment *other) const
unsigned int getEmbeddedNodeIndex(EFANode *node) const
bool isFaceInterior(unsigned int face_id) const
const dof_id_type n_nodes
std::vector< EFANode * > _nodes
std::vector< EFANode * > getCommonNodes(const EFAElement3D *other_elem) const
EFAFace * getFace(unsigned int face_id) const
unsigned int numNodes() const
virtual bool willCrackTipExtend(std::vector< unsigned int > &split_neighbors) const
virtual void getNeighborEdgeIndex(const EFAElement3D *neighbor_elem, unsigned int face_id, unsigned int edge_id, unsigned int &neigh_face_id, unsigned int &neigh_edge_id) const
void setParent(EFAElement *parent)
std::vector< EFAFragment3D * > _fragments
void findGeneralNeighbors(const std::map< EFANode *, std::set< EFAElement *>> &InverseConnectivity)
virtual void removePhantomEmbeddedNode()
bool containsEdge(const EFAEdge &other) const
std::set< EFANode * > getTipEmbeddedNodes() const
double linearTetShape3D(unsigned int node_id, std::vector< double > &xi_3d)
virtual bool isPartial() const
EFAEdge * getEdge(unsigned int edge_id) const
virtual void switchNode(EFANode *new_node, EFANode *old_node, bool descend_to_parent)=0
virtual void switchEmbeddedNode(EFANode *new_node, EFANode *old_node)
std::vector< unsigned int > getTipFaceIDs() const
EFANode * getNode(unsigned int node_id) const
unsigned int numCommonElems(std::set< T > &v1, std::set< T > &v2)
bool equivalent(const EFAEdge &other) const
double linearHexShape3D(unsigned int node_id, std::vector< double > &xi_3d)
virtual bool shouldDuplicateCrackTipSplitElement(const std::set< EFAElement *> &CrackTipElements)
unsigned int _num_vertices
EFANode * getInteriorFaceNode(unsigned int i) const
EFAElement * getGeneralNeighbor(unsigned int index) const
virtual void neighborSanityCheck() const
bool getFragmentFaceID(unsigned int elem_face_id, unsigned int &frag_face_id) const
bool fragmentFaceAlreadyCut(unsigned int ElemFaceID) const
virtual void restoreFragment(const EFAElement *const from_elem)
virtual void getMasterInfo(EFANode *node, std::vector< EFANode *> &master_nodes, std::vector< double > &master_weights) const
void removeEmbeddedNode(EFANode *emb_node, bool remove_for_neighbor)
void mergeNodes(EFANode *&childNode, EFANode *&childOfNeighborNode, EFAElement *childOfNeighborElem, std::map< unsigned int, EFANode *> &PermanentNodes, std::map< unsigned int, EFANode *> &TempNodes)
bool hasIntersectionAtPosition(double position, EFANode *from_node) const
bool isPhysicalEdgeCut(unsigned int ElemFaceID, unsigned int ElemFaceEdgeID, double position) const
std::vector< unsigned int > getCommonFaceID(const EFAElement3D *other_elem) const
virtual void createChild(const std::set< EFAElement *> &CrackTipElements, std::map< unsigned int, EFAElement *> &Elements, std::map< unsigned int, EFAElement *> &newChildElements, std::vector< EFAElement *> &ChildElements, std::vector< EFAElement *> &ParentElements, std::map< unsigned int, EFANode *> &TempNodes)
virtual void setupNeighbors(std::map< EFANode *, std::set< EFAElement *>> &InverseConnectivityMap)
std::vector< T > getCommonElems(std::set< T > &v1, std::set< T > &v2)
unsigned int getNeighborFaceInteriorNodeID(unsigned int face_id, unsigned int node_id, EFAElement3D *neighbor_elem) const
IntRange< T > make_range(T beg, T end)
EFANode * createLocalNodeFromGlobalNode(const EFANode *global_node) const
virtual void updateFragmentNode()
bool overlaysElement(const EFAElement3D *other_elem) const
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
void normalizePoint(Point &p)
unsigned int numNodes() const
unsigned int numEmbeddedNodes() const
void addIntersection(double position, EFANode *embedded_node_tmp, EFANode *from_node)
void setLocalCoordinates()
EFAElement3D(unsigned int eid, unsigned int n_nodes, unsigned int n_faces)
unsigned int adjacentCommonEdge(const EFAFace *other_face) const
EFAFace * getFragmentFace(unsigned int frag_id, unsigned int face_id) const
std::vector< EFAElement * > _general_neighbors
std::vector< EFAElement * > _children
void addFragFaceEdgeCut(unsigned int frag_face_id, unsigned int frag_edge_id, double position, std::map< unsigned int, EFANode *> &EmbeddedNodes, bool add_to_neighbor, bool add_to_adjacent)
bool _crack_tip_split_element
EFANode * getNode(unsigned int index) const
virtual void fragmentSanityCheck(unsigned int n_old_frag_faces, unsigned int n_old_frag_cuts) const
std::vector< unsigned int > _crack_tip_neighbors
static const std::string k
unsigned int numFaceNeighbors(unsigned int face_id) const
virtual unsigned int getNumCuts() const
auto index_range(const T &sizable)
std::vector< EFAVolumeNode * > _interior_nodes
EFAVolumeNode * getInteriorNode(unsigned int interior_node_id) const
static const std::string center
void setNode(unsigned int node_id, EFANode *node)
void addFaceEdgeCut(unsigned int face_id, unsigned int edge_id, double position, EFANode *embedded_node, std::map< unsigned int, EFANode *> &EmbeddedNodes, bool add_to_neighbor, bool add_to_adjacent)
virtual void clearNeighbors()
std::set< EFANode * > getFaceNodes(unsigned int face_id) const
virtual bool isFinalCut() const
bool faceContainsTip(unsigned int face_id) const
std::vector< EFANode * > _local_nodes
unsigned int numFaces() const
std::set< EFANode * > getPhantomNodeOnFace(unsigned int face_id) const
N_CATEGORY category() const
bool fragmentHasTipFaces() const
virtual bool isCrackTipElement() const