Skip to content

Commit 886eb49

Browse files
committed
Update boundary_info during poly2tri refinement
This fixes a number of bugs that trigger for me in complicated refinement cases.
1 parent ff71886 commit 886eb49

1 file changed

Lines changed: 46 additions & 13 deletions

File tree

src/mesh/poly2tri_triangulator.C

Lines changed: 46 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -806,7 +806,8 @@ bool Poly2TriTriangulator::insert_refinement_points()
806806
// cavity will be a source of one new triangle.
807807

808808
// Keep maps for doing neighbor pointer assignment
809-
std::map<Node *, Elem *> neighbors_CCW, neighbors_CW;
809+
std::map<Node *, std::pair<Elem *, boundary_id_type>>
810+
neighbors_CCW, neighbors_CW;
810811

811812
for (Elem * old_elem : cavity)
812813
{
@@ -818,23 +819,35 @@ bool Poly2TriTriangulator::insert_refinement_points()
818819
Node * node_CW = old_elem->node_ptr(s),
819820
* node_CCW = old_elem->node_ptr((s+1)%3);
820821

821-
auto set_neighbors = [&neighbors_CW, &neighbors_CCW,
822-
&node_CW, &node_CCW](Elem * new_neigh)
822+
auto set_neighbors =
823+
[&neighbors_CW, &neighbors_CCW, &node_CW,
824+
&node_CCW, &boundary_info]
825+
(Elem * new_neigh, boundary_id_type bcid)
823826
{
824827
// Set clockwise neighbor and vice-versa if possible
825828
auto CW_it = neighbors_CW.find(node_CW);
826829
if (CW_it == neighbors_CW.end())
827830
{
828831
libmesh_assert(!neighbors_CCW.count(node_CW));
829-
neighbors_CCW[node_CW] = new_neigh;
832+
neighbors_CCW[node_CW] = std::make_pair(new_neigh, bcid);
830833
}
831834
else
832835
{
833-
Elem * neigh_CW = CW_it->second;
836+
Elem * neigh_CW = CW_it->second.first;
834837
if (new_neigh)
835-
new_neigh->set_neighbor(0, neigh_CW);
838+
{
839+
new_neigh->set_neighbor(0, neigh_CW);
840+
boundary_id_type bcid_CW = CW_it->second.second;
841+
if (bcid_CW != BoundaryInfo::invalid_id)
842+
boundary_info.add_side(new_neigh, 0, bcid_CW);
843+
844+
}
836845
if (neigh_CW)
837-
neigh_CW->set_neighbor(2, new_neigh);
846+
{
847+
neigh_CW->set_neighbor(2, new_neigh);
848+
if (bcid != BoundaryInfo::invalid_id)
849+
boundary_info.add_side(neigh_CW, 2, bcid);
850+
}
838851
neighbors_CW.erase(CW_it);
839852
}
840853

@@ -843,15 +856,24 @@ bool Poly2TriTriangulator::insert_refinement_points()
843856
if (CCW_it == neighbors_CCW.end())
844857
{
845858
libmesh_assert(!neighbors_CW.count(node_CCW));
846-
neighbors_CW[node_CCW] = new_neigh;
859+
neighbors_CW[node_CCW] = std::make_pair(new_neigh, bcid);
847860
}
848861
else
849862
{
850-
Elem * neigh_CCW = CCW_it->second;
863+
Elem * neigh_CCW = CCW_it->second.first;
851864
if (new_neigh)
852-
new_neigh->set_neighbor(2, neigh_CCW);
865+
{
866+
boundary_id_type bcid_CCW = CCW_it->second.second;
867+
new_neigh->set_neighbor(2, neigh_CCW);
868+
if (bcid_CCW != BoundaryInfo::invalid_id)
869+
boundary_info.add_side(new_neigh, 2, bcid_CCW);
870+
}
853871
if (neigh_CCW)
854-
neigh_CCW->set_neighbor(0, new_neigh);
872+
{
873+
neigh_CCW->set_neighbor(0, new_neigh);
874+
if (bcid != BoundaryInfo::invalid_id)
875+
boundary_info.add_side(neigh_CCW, 0, bcid);
876+
}
855877
neighbors_CCW.erase(CCW_it);
856878
}
857879
};
@@ -862,7 +884,10 @@ bool Poly2TriTriangulator::insert_refinement_points()
862884
if (old_elem == cavity_elem &&
863885
s == side)
864886
{
865-
set_neighbors(nullptr);
887+
std::vector<boundary_id_type> bcids;
888+
boundary_info.boundary_ids(old_elem, s, bcids);
889+
libmesh_assert_equal_to(bcids.size(), 1);
890+
set_neighbors(nullptr, bcids[0]);
866891
continue;
867892
}
868893

@@ -879,9 +904,15 @@ bool Poly2TriTriangulator::insert_refinement_points()
879904
neigh->which_neighbor_am_i(old_elem);
880905
neigh->set_neighbor(neigh_s, new_elem.get());
881906
}
907+
else
908+
{
909+
std::vector<boundary_id_type> bcids;
910+
boundary_info.boundary_ids(old_elem, s, bcids);
911+
boundary_info.add_side(new_elem.get(), 1, bcids);
912+
}
882913

883914
// Set in-cavity neighbors' neighbor pointers
884-
set_neighbors(new_elem.get());
915+
set_neighbors(new_elem.get(), BoundaryInfo::invalid_id);
885916

886917
// C++ allows function argument evaluation in any
887918
// order, but we need get() to precede move
@@ -890,6 +921,8 @@ bool Poly2TriTriangulator::insert_refinement_points()
890921
}
891922
}
892923

924+
boundary_info.remove(old_elem);
925+
893926
auto it = new_elems.find(old_elem);
894927
if (it == new_elems.end())
895928
mesh.delete_elem(old_elem);

0 commit comments

Comments
 (0)