mesh_graph/ops/cleanup/

vertex_neighborhood.rs

use hashbrown::HashSet;
use itertools::Itertools;
use tracing::{error, instrument};

use crate::{FaceId, HalfedgeId, MeshGraph, VertexId, error_none};

#[derive(Default)]
pub struct VertexNeighborhoodCleanup {
    pub removed_vertices: Vec<VertexId>,
    pub removed_halfedges: Vec<HalfedgeId>,
    pub removed_faces: Vec<FaceId>,

    pub added_vertices: Vec<VertexId>,
}

#[derive(Default)]
struct VertexNeighborhoodCleanupStep {
    pub removed_vertices: Vec<VertexId>,
    pub removed_halfedges: Vec<HalfedgeId>,
    pub removed_faces: Vec<FaceId>,

    pub added_duplicated_vertices: Vec<VertexId>,
    pub touched_vertices: Vec<VertexId>,
}

impl MeshGraph {
    /// Ensure that the neighborhood of a vertex is manifold.
    ///
    /// It removes flaps (two neighboring coincident triangles) and splits the given vertex
    /// in two if there are non-neighboring degenerate triangles or edges.
    ///
    /// See [Freestyle: Sculpting meshes with self-adaptive topology DOI 10.1016/j.cag.2011.03.033](https://inria.hal.science/inria-00606516v1/document)
    /// Chapters 3.2 and 5.1
    #[instrument(skip(self))]
    pub fn make_vertex_neighborhood_manifold(
        &mut self,
        vertex_id: VertexId,
    ) -> VertexNeighborhoodCleanup {
        #[cfg(feature = "rerun")]
        self.log_vert_rerun("make_neigh_manifold", vertex_id);

        let mut result = VertexNeighborhoodCleanup::default();

        let mut vertices = vec![vertex_id];

        let vertex_count = self.vertices.len();
        let mut iterations = 0;

        while let Some(&v_id) = vertices.first()
            && iterations < vertex_count
        {
            if !self.vertices.contains_key(v_id) {
                vertices.swap_remove(0);
                continue;
            }

            let VertexNeighborhoodCleanupStep {
                added_duplicated_vertices,
                touched_vertices,
                removed_vertices,
                removed_halfedges,
                removed_faces,
            } = self.make_vertex_neighborhood_manifold_step(v_id);

            if removed_vertices.is_empty()
                && removed_halfedges.is_empty()
                && removed_faces.is_empty()
                && added_duplicated_vertices.is_empty()
                && touched_vertices.is_empty()
            {
                vertices.swap_remove(0);
            }

            vertices.extend(touched_vertices);
            for added_vertex in added_duplicated_vertices {
                vertices.push(added_vertex);

                result.added_vertices.push(added_vertex);
            }

            result.removed_vertices.extend(removed_vertices);
            result.removed_halfedges.extend(removed_halfedges);
            result.removed_faces.extend(removed_faces);

            iterations += 1;
        }

        for &cancelled_v_id in
            HashSet::<VertexId>::from_iter(result.removed_vertices.iter().copied())
                .intersection(&HashSet::from_iter(result.added_vertices.iter().copied()))
        {
            result.added_vertices.retain(|&v_id| v_id != cancelled_v_id);
            result
                .removed_vertices
                .retain(|&v_id| v_id != cancelled_v_id);
        }

        result
    }

    fn make_vertex_neighborhood_manifold_step(
        &mut self,
        vertex_id: VertexId,
    ) -> VertexNeighborhoodCleanupStep {
        let mut removed_vertices = Vec::new();
        let mut removed_halfedges = Vec::new();
        let mut removed_faces = Vec::new();

        let added_duplicated_vertices = self
            .make_vertex_neighborhood_manifold_inner(
                vertex_id,
                &mut removed_vertices,
                &mut removed_halfedges,
                &mut removed_faces,
            )
            .unwrap_or_default();

        VertexNeighborhoodCleanupStep {
            added_duplicated_vertices,
            touched_vertices: vec![],
            removed_vertices,
            removed_halfedges,
            removed_faces,
        }
    }

    pub fn make_vertex_neighborhood_manifold_inner(
        &mut self,
        vertex_id: VertexId,
        removed_vertices: &mut Vec<VertexId>,
        removed_halfedges: &mut Vec<HalfedgeId>,
        removed_faces: &mut Vec<FaceId>,
    ) -> Option<Vec<VertexId>> {
        self.vertices.get(vertex_id)?;

        if self.remove_single_face(
            vertex_id,
            removed_vertices,
            removed_halfedges,
            removed_faces,
        )? {
            return Some(vec![]);
        }

        if self.remove_neighboring_flaps(
            vertex_id,
            removed_vertices,
            removed_halfedges,
            removed_faces,
        )? {
            return Some(vec![]);
        }

        if let Some(new_vertices) = self.split_disconnected_neighborhoods(vertex_id)
            && !new_vertices.is_empty()
        {
            return Some(new_vertices);
        }

        if let Some(inserted_duplicated_vertex) = self.remove_degenerate_faces(
            vertex_id,
            removed_vertices,
            removed_halfedges,
            removed_faces,
        ) {
            return Some(vec![inserted_duplicated_vertex]);
        }

        self.remove_degenerate_edges(vertex_id).map(|v| vec![v])
    }

    fn split_disconnected_neighborhoods(&mut self, vertex_id: VertexId) -> Option<Vec<VertexId>> {
        let mut new_vertices = Vec::new();

        let mut outgoing_halfedges = HashSet::<HalfedgeId>::from_iter(
            self.outgoing_halfedges
                .get(vertex_id)
                .or_else(error_none!("Outgoing halfedges not found"))?
                .iter()
                .copied(),
        );

        while let Some(&start_he_id) = outgoing_halfedges.iter().next() {
            let mut current_he_id = start_he_id;

            let len = outgoing_halfedges.len();
            let mut current_outgoing_halfedges = Vec::with_capacity(len);

            loop {
                if !outgoing_halfedges.remove(&current_he_id) {
                    return None;
                }
                current_outgoing_halfedges.push(current_he_id);
                if current_outgoing_halfedges.len() == len {
                    return None;
                }

                let cur_he = self
                    .halfedges
                    .get(current_he_id)
                    .or_else(error_none!("Halfedge not found"))?;

                let Some(next_he_id) = cur_he.cw_rotated_neighbour(self) else {
                    // The fan runs into a boundary, so rotating clockwise can't
                    // reach the rest of it. Whether the walk started in the
                    // middle of the fan or at its end decides how much is left,
                    // so collect the remainder by rotating the other way around
                    // the vertex, back from `start_he_id`. Otherwise a perfectly
                    // connected boundary neighborhood looks disconnected and the
                    // vertex gets split for no reason.
                    let mut back_he_id = start_he_id;
                    while let Some(prev_he_id) = self
                        .halfedges
                        .get(back_he_id)
                        .or_else(error_none!("Halfedge not found"))?
                        .ccw_rotated_neighbour(self)
                    {
                        if !outgoing_halfedges.remove(&prev_he_id) {
                            break;
                        }
                        current_outgoing_halfedges.push(prev_he_id);
                        back_he_id = prev_he_id;
                    }

                    break;
                };

                if next_he_id == start_he_id {
                    break;
                }

                current_he_id = next_he_id;
            }

            if !outgoing_halfedges.is_empty() {
                let new_vertex_id = self
                    .duplicate_vertex_and_assign_halfedges(vertex_id, current_outgoing_halfedges)?;

                new_vertices.push(new_vertex_id);

                self.vertices
                    .get_mut(vertex_id)
                    .or_else(error_none!("Vertex not found"))?
                    .outgoing_halfedge = Some(*outgoing_halfedges.iter().next().unwrap()); // check above that outgoing_halfedges is not empty

                self.outgoing_halfedges
                    .insert(vertex_id, Vec::from_iter(outgoing_halfedges));

                return Some(new_vertices);
            }
        }

        Some(new_vertices)
    }

    fn duplicate_vertex_and_assign_halfedges(
        &mut self,
        vertex_id: VertexId,
        outgoing_halfedges: Vec<HalfedgeId>,
    ) -> Option<VertexId> {
        let new_vert_id = self.add_vertex(
            *self
                .positions
                .get(vertex_id)
                .or_else(error_none!("Position not found"))?
                + glam::Vec3::new(0.1, 0.0, 0.0),
        );

        tracing::debug!("Duplicated {vertex_id:?}: {new_vert_id:?}");

        for &he_id in &outgoing_halfedges {
            let he = self
                .halfedges
                .get(he_id)
                .or_else(error_none!("Halfedge not found"))?;

            let twin_id = he.twin.or_else(error_none!("Twin not found"))?;

            self.halfedges
                .get_mut(twin_id)
                .or_else(error_none!("Twin not found"))?
                .end_vertex = new_vert_id;
        }

        // just created above
        self.vertices[new_vert_id].outgoing_halfedge = Some(outgoing_halfedges[0]);

        self.outgoing_halfedges
            .insert(new_vert_id, outgoing_halfedges);

        Some(new_vert_id)
    }

    fn remove_degenerate_edges(&mut self, vertex_id: VertexId) -> Option<VertexId> {
        let halfedges = self
            .vertices
            .get(vertex_id)
            .or_else(error_none!("Vertex not found"))?
            .outgoing_halfedges(self)
            .collect_vec();

        for [he_id1, he_id2] in halfedges.into_iter().array_combinations() {
            if self.halfedges_share_all_vertices(he_id1, he_id2) {
                let he1 = self.halfedges[he_id1]; // checked in `halfedges_share_all_vertices`
                let twin_id1 = he1.twin.or_else(error_none!("Twin not found"))?;
                let face_id1 = self
                    .halfedges
                    .get(twin_id1)
                    .or_else(error_none!("Halfedge not found"))?
                    .face
                    .or_else(error_none!("Face not found"))?;
                let face_id2 = self.halfedges[he_id2] // checked in `halfedges_share_all_vertices`
                    .face
                    .or_else(error_none!("Face not found"))?;

                // checked above
                let coincident_face_ids = self.vertices[vertex_id].faces(self).collect_vec();
                let mut start_idx = 0;

                for (idx, face_id) in coincident_face_ids.iter().enumerate() {
                    if *face_id == face_id1 {
                        start_idx = idx;
                        break;
                    }
                }

                let mut end_idx = start_idx;
                let mut side_one = vec![];

                for _ in 0..coincident_face_ids.len() {
                    let face_id = coincident_face_ids[end_idx];

                    side_one.push(face_id);

                    end_idx += 1;
                    end_idx %= coincident_face_ids.len();

                    if face_id == face_id2 {
                        break;
                    }
                }

                let mut side_two = vec![];

                for _ in 0..coincident_face_ids.len() {
                    if end_idx == start_idx {
                        break;
                    }

                    let face_id = coincident_face_ids[end_idx];

                    side_two.push(face_id);

                    end_idx += 1;
                    end_idx %= coincident_face_ids.len();
                }

                return self.split_regions_at_edge(vertex_id, he1.end_vertex, side_one, side_two);
            }
        }

        None
    }

    #[instrument(skip(self))]
    fn split_regions_at_edge(
        &mut self,
        vertex_id: VertexId,
        other_vertex_id: VertexId,
        side_one: Vec<FaceId>,
        side_two: Vec<FaceId>,
    ) -> Option<VertexId> {
        if side_one.len() < 2 || side_two.len() < 2 {
            error!("Not enough halfedges to split");
            return None;
        }

        #[cfg(feature = "rerun")]
        {
            self.log_vert_rerun("split_regions_at_edge", vertex_id);
            self.log_faces_rerun("split_regions_at_edge/side_one", &side_one);
            self.log_faces_rerun("split_regions_at_edge/side_two", &side_two);
        }

        let new_vertex_id = self.add_vertex(
            *self
                .positions
                .get(vertex_id)
                .or_else(error_none!("Vertex position not found"))?,
        );

        // Updating vertices
        for face_id in &side_two {
            let face = self
                .faces
                .get(*face_id)
                .or_else(error_none!("Face not found"))?;

            for he_id in face.halfedges(self).collect_vec() {
                // already checked in iterator that this he exists
                let he = &mut self.halfedges[he_id];

                if he.end_vertex == vertex_id {
                    he.end_vertex = new_vertex_id;
                }
            }
        }

        self.weld_faces_at(
            vertex_id,
            other_vertex_id,
            side_one[side_one.len() - 1],
            side_one[0],
        );
        self.weld_faces_at(
            new_vertex_id,
            other_vertex_id,
            side_two[0],
            side_two[side_two.len() - 1],
        );

        // Ground-truth rebuild: a seed-walk rebuild could miss halfedges behind a
        // temporarily detached fan (the re-points above already partitioned the
        // star), which would silently drop live outgoing entries.
        #[cfg(feature = "instrumentation")]
        crate::record_op_trace!(
            "split_regions_at_edge({vertex_id:?}, {other_vertex_id:?}): new {new_vertex_id:?}"
        );
        self.rebuild_vertex_outgoing_list(vertex_id);
        self.rebuild_vertex_outgoing_list(new_vertex_id);

        Some(new_vertex_id)
    }

    #[instrument(skip_all)]
    pub fn remove_degenerate_faces(
        &mut self,
        vertex_id: VertexId,
        removed_vertices: &mut Vec<VertexId>,
        removed_halfedges: &mut Vec<HalfedgeId>,
        removed_faces: &mut Vec<FaceId>,
    ) -> Option<VertexId> {
        let faces = self
            .vertices
            .get(vertex_id)
            .or_else(error_none!("Vertex not found"))?
            .faces(self)
            .collect_vec();

        for [face_id1, face_id2] in faces.into_iter().array_combinations() {
            if self.faces_share_all_vertices(face_id1, face_id2) {
                // checked above
                let coincident_face_ids = self.vertices[vertex_id].faces(self).collect_vec();
                let mut start_idx = 0;

                for (idx, face_id) in coincident_face_ids.iter().enumerate() {
                    if *face_id == face_id1 {
                        start_idx = (idx + 1) % coincident_face_ids.len();
                        break;
                    }
                }

                let mut end_idx = start_idx;
                let mut side_one = vec![];

                for _ in 0..coincident_face_ids.len() {
                    let face_id = coincident_face_ids[end_idx];

                    end_idx += 1;
                    end_idx %= coincident_face_ids.len();

                    if face_id == face_id2 {
                        break;
                    }

                    side_one.push(face_id);
                }

                let mut side_two = vec![];

                for _ in 0..coincident_face_ids.len() {
                    if end_idx == start_idx {
                        break;
                    }

                    let face_id = coincident_face_ids[end_idx];

                    side_two.push(face_id);

                    end_idx += 1;
                    end_idx %= coincident_face_ids.len();
                }

                side_two.pop();

                let new_vertex_id =
                    self.split_regions_at_vertex(vertex_id, side_one, side_two, removed_halfedges);

                let (del_v_ids, del_he_ids) = self.remove_face(face_id1);
                removed_vertices.extend(del_v_ids);
                removed_halfedges.extend(del_he_ids);
                removed_faces.push(face_id1);

                let (del_v_ids, del_he_ids) = self.remove_face(face_id2);
                removed_vertices.extend(del_v_ids);
                removed_halfedges.extend(del_he_ids);
                removed_faces.push(face_id2);

                return new_vertex_id;
            }
        }

        None
    }

    #[instrument(skip(self))]
    fn split_regions_at_vertex(
        &mut self,
        vertex_id: VertexId,
        side_one: Vec<FaceId>,
        side_two: Vec<FaceId>,
        removed_halfedges: &mut Vec<HalfedgeId>,
    ) -> Option<VertexId> {
        if side_one.len() < 2 || side_two.len() < 2 {
            error!("Not enough halfedges to split");
            return None;
        }

        #[cfg(feature = "rerun")]
        {
            self.log_vert_rerun("split_regions_at_vertex", vertex_id);
            self.log_faces_rerun("split_regions_at_vertex/side_one", &side_one);
            self.log_faces_rerun("split_regions_at_vertex/side_two", &side_two);
        }

        let new_vertex_id = self.add_vertex(
            *self
                .positions
                .get(vertex_id)
                .or_else(error_none!("Vertex position not found"))?,
        );

        // Updating start vertex
        for face_id in &side_two {
            let face = self
                .faces
                .get(*face_id)
                .or_else(error_none!("Face not found"))?;

            for he_id in face.halfedges(self).collect_vec() {
                // already checked in iterator that this he exists
                let he = &mut self.halfedges[he_id];
                if he.end_vertex == vertex_id {
                    he.end_vertex = new_vertex_id;
                }
            }
        }

        self.weld_faces(vertex_id, side_one[side_one.len() - 1], side_one[0]);
        self.weld_faces(new_vertex_id, side_two[0], side_two[side_two.len() - 1]);

        // Ground-truth rebuild: see `split_regions_at_edge`.
        #[cfg(feature = "instrumentation")]
        crate::record_op_trace!("split_regions_at_vertex({vertex_id:?}): new {new_vertex_id:?}");
        self.rebuild_vertex_outgoing_list(vertex_id);
        self.rebuild_vertex_outgoing_list(new_vertex_id);

        Some(new_vertex_id)
    }

    /// This finds the common edge between the two faces and welds them together by connecting
    /// the two `twin` relationships. The old twins are re-paired with each other so
    /// every halfedge stays paired and the twin relationship stays symmetric (see
    /// [`Self::weld_faces_at`]).
    ///
    /// `start_vertex_id` is shared by `face_id1` and `face_id2`.
    fn weld_faces(
        &mut self,
        start_vertex_id: VertexId,
        face_id1: FaceId,
        face_id2: FaceId,
    ) -> Option<()> {
        // #[cfg(feature = "rerun")]
        // {
        //     self.log_face_rerun("face1", face_id1);
        //     self.log_face_rerun("face2", face_id2);
        // }

        let face1 = self
            .faces
            .get(face_id1)
            .or_else(error_none!("Face not found"))?;
        let face2 = self
            .faces
            .get(face_id2)
            .or_else(error_none!("Face not found"))?;

        let face1_vertices = face1.vertices(self).collect::<HashSet<_>>();
        let face2_vertices = face2.vertices(self).collect::<HashSet<_>>();

        let mut other_common_vertex_id = None;

        for &v_id in face1_vertices.intersection(&face2_vertices) {
            if v_id != start_vertex_id {
                other_common_vertex_id = Some(v_id);
                break;
            }
        }

        let other_common_vertex_id =
            other_common_vertex_id.or_else(error_none!("No other common vertex found"))?;

        self.weld_faces_at(start_vertex_id, other_common_vertex_id, face_id1, face_id2)
    }

    #[instrument(skip(self))]
    fn weld_faces_at(
        &mut self,
        start_vertex_id: VertexId,
        other_common_vertex_id: VertexId,
        face_id1: FaceId,
        face_id2: FaceId,
    ) -> Option<()> {
        let face1 = self
            .faces
            .get(face_id1)
            .or_else(error_none!("Face not found"))?;
        let face2 = self
            .faces
            .get(face_id2)
            .or_else(error_none!("Face not found"))?;

        let he1_id = face1
            .halfedge_between(start_vertex_id, other_common_vertex_id, self)
            .or_else(error_none!("Halfedge between vertices not found"))?;

        let he2_id = face2
            .halfedge_between(start_vertex_id, other_common_vertex_id, self)
            .or_else(error_none!("Halfedge between vertices not found"))?;

        // The face chains can be stale after flap removals, so verify the halfedges are
        // still present before re-pairing them (a dangling twin write corrupts the mesh).
        if !self.halfedges.contains_key(he1_id) || !self.halfedges.contains_key(he2_id) {
            error!("Halfedge not found in weld_faces_at");
            return None;
        }
        // Two coincident faces on the same chain share the same halfedge(s): re-pairing
        // would write a self-twin (`he1.twin = he1`), which later makes the collapse's
        // `remove_halfedge_face` re-pair fail and strands live face members twinless.
        // Leave the pairing untouched; coincident faces are handled by the flap cleanup.
        if he1_id == he2_id {
            error!("weld_faces_at: faces share the same halfedge {he1_id:?}");
            return None;
        }
        // The old twins' back-references must not be left dangling: if `he1_id` was not
        // already `he2_id`'s twin, the weld re-pairs the two faces' edge halves while
        // the old partners (`he1_id`'s and `he2_id`'s previous twins) would keep
        // pointing at the re-paired halfedges — one-sided references that later
        // removals cannot see from the removed halfedge's own twin field and that no
        // cleanup pass repairs anymore. When both old partners are clean symmetric
        // partners they are re-paired with each other, which keeps every halfedge
        // paired (no twinless survivors, no one-sided twins); a partner that cannot
        // be re-paired is severed instead.
        let t1 = self.halfedges.get(he1_id).and_then(|h| h.twin);
        let t2 = self.halfedges.get(he2_id).and_then(|h| h.twin);
        let t1_swappable = t1.is_some_and(|t1| {
            t1 != he2_id
                && self
                    .halfedges
                    .get(t1)
                    .is_some_and(|t| t.twin == Some(he1_id))
        });
        let t2_swappable = t2.is_some_and(|t2| {
            t2 != he1_id
                && self
                    .halfedges
                    .get(t2)
                    .is_some_and(|t| t.twin == Some(he2_id))
        });
        match (t1_swappable, t2_swappable) {
            (true, true) => {
                // Both old partners are clean: re-pair them with each other (the two
                // duplicate edges swap partners), keeping all four halfedges paired.
                let t1 = t1.unwrap();
                let t2 = t2.unwrap();
                if t1 == t2 {
                    // Both sides are the same halfedge (identical faces); pairing it
                    // with itself is a self-twin. Leave it as-is.
                    error!("weld_faces_at: old partners are identical {t1:?}");
                } else {
                    self.halfedges.get_mut(t1).unwrap().twin = Some(t2);
                    self.halfedges.get_mut(t2).unwrap().twin = Some(t1);
                }
            }
            (true, false) => {
                // `he2`'s old partner is gone/unpaired: `t1` loses its partner when
                // `he1` is re-paired below. Give it a fresh boundary partner instead
                // of leaving it twinless.
                let t1 = t1.unwrap();
                // `t1` is `he1`'s twin, so `t1.start == he1.end_vertex`.
                let t1_start = self
                    .halfedges
                    .get(he1_id)
                    .map(|h| h.end_vertex)
                    .or_else(error_none!("Weld sever start vertex not found"))?;
                if self.pair_with_fresh_boundary_half(t1, t1_start).is_none() {
                    error!("weld_faces_at: could not re-pair severed partner {t1:?}");
                }
            }
            (false, true) => {
                // `t2` is `he2`'s twin, so `t2.start == he2.end_vertex`.
                let t2 = t2.unwrap();
                let t2_start = self
                    .halfedges
                    .get(he2_id)
                    .map(|h| h.end_vertex)
                    .or_else(error_none!("Weld sever start vertex not found"))?;
                if self.pair_with_fresh_boundary_half(t2, t2_start).is_none() {
                    error!("weld_faces_at: could not re-pair severed partner {t2:?}");
                }
            }
            _ => {}
        }
        self.halfedges[he1_id].twin = Some(he2_id);
        self.halfedges[he2_id].twin = Some(he1_id);

        let (start_out_he_id, other_out_he_id) =
            if self.halfedges[he1_id].end_vertex == start_vertex_id {
                (he2_id, he1_id)
            } else {
                (he1_id, he2_id)
            };

        self.vertices
            .get_mut(start_vertex_id)
            .or_else(error_none!("Start vertex not found"))?
            .outgoing_halfedge = Some(start_out_he_id);
        self.vertices
            .get_mut(other_common_vertex_id)
            .or_else(error_none!("Other vertex not found"))?
            .outgoing_halfedge = Some(other_out_he_id);

        Some(())
    }

    pub(crate) fn remove_single_face(
        &mut self,
        vertex_id: VertexId,
        removed_vertices: &mut Vec<VertexId>,
        removed_halfedges: &mut Vec<HalfedgeId>,
        removed_faces: &mut Vec<FaceId>,
    ) -> Option<bool> {
        let face_ids = self
            .outgoing_halfedges
            .get(vertex_id)
            .or_else(error_none!(
                "Outgoing halfedges not found for vertex {vertex_id:?}"
            ))?
            .iter()
            .filter_map(|&he_id| {
                self.halfedges
                    .get(he_id)
                    .or_else(error_none!("Halfedge not found"))
                    .and_then(|he| he.face)
            })
            .collect_vec();

        if face_ids.len() == 1 {
            let face_id = face_ids[0];

            let (v_ids, he_ids) = self.remove_face(face_id);

            removed_vertices.extend(v_ids);
            removed_halfedges.extend(he_ids);
            removed_faces.push(face_id);

            Some(true)
        } else {
            Some(false)
        }
    }

    /// Removes all neighboring flaps (triangles that share all the same vertices) connected to the given vertex.
    pub fn remove_neighboring_flaps(
        &mut self,
        vertex_id: VertexId,
        removed_vertices: &mut Vec<VertexId>,
        removed_halfedges: &mut Vec<HalfedgeId>,
        removed_faces: &mut Vec<FaceId>,
    ) -> Option<bool> {
        let faces = self
            .vertices
            .get(vertex_id)
            .or_else(error_none!("Vertex not found"))?
            .faces(self)
            .collect_vec();

        if faces.len() < 2 {
            return Some(false);
        }

        let mut face_tuples = faces.into_iter().circular_array_windows().collect_vec();

        while let Some([face_id1, face_id2]) = face_tuples.pop() {
            if self.faces_share_all_vertices(face_id1, face_id2) {
                #[cfg(feature = "rerun")]
                {
                    self.log_vert_rerun("flap", vertex_id);
                    self.log_face_rerun("flap1", face_id1);
                    self.log_face_rerun("flap2", face_id2);
                }

                let mut halfedges_of_faces = HashSet::<HalfedgeId>::from_iter(
                    self.halfedges.iter().filter_map(|(he_id, he)| {
                        if he.face == Some(face_id1) || he.face == Some(face_id2) {
                            Some(he_id)
                        } else {
                            None
                        }
                    }),
                );

                let (vs, hes) = self.remove_face(face_id1);
                for he_id in &hes {
                    halfedges_of_faces.remove(he_id);
                }
                removed_vertices.extend(vs);
                removed_halfedges.extend(hes);
                removed_faces.push(face_id1);

                let (vs, hes) = self.remove_face(face_id2);
                for he_id in &hes {
                    halfedges_of_faces.remove(he_id);
                }
                removed_vertices.extend(vs);
                removed_halfedges.extend(hes);
                removed_faces.push(face_id2);

                // TODO : Handle the case when there is only one halfedge?
                if halfedges_of_faces.len() >= 2 {
                    for [he_id1, he_id2] in halfedges_of_faces.iter().copied().array_combinations()
                    {
                        let Some(he1) = self.halfedges.get(he_id1) else {
                            // might have been deleted by prior iteration
                            continue;
                        };
                        let twin_id1 = he1.twin.or_else(error_none!("Twin 1 missing"))?;

                        let Some(he2) = self.halfedges.get(he_id2) else {
                            continue;
                        };
                        let twin_id2 = he2.twin.or_else(error_none!("Twin 2 missing"))?;

                        let start_v_id1 = he1
                            .start_vertex(self)
                            .or_else(error_none!("Start vertex 1 missing"))?;
                        let start_v_id2 = he2
                            .start_vertex(self)
                            .or_else(error_none!("Start vertex 2 missing"))?;

                        if he1.end_vertex != start_v_id2 || he2.end_vertex != start_v_id1 {
                            continue;
                        }

                        // the twin edges of the neighboring faces of the deleted faces are still there
                        // we need to remove them and re-connect (twin) their twin edges

                        let he1_end_v = he1.end_vertex;
                        let he2_end_v = he2.end_vertex;

                        self.remove_only_halfedge(he_id1);
                        self.remove_only_halfedge(he_id2);
                        removed_halfedges.push(he_id1);
                        removed_halfedges.push(he_id2);

                        #[cfg(feature = "instrumentation")]
                        crate::record_op_trace!(
                            "flap removal: removed {he_id1:?}+{he_id2:?} (edge {start_v_id1:?}-{start_v_id2:?}); twins {twin_id1:?}/{twin_id2:?}"
                        );

                        // `remove_only_halfedge` above cleared the twins' back-pointers, so
                        // `twin_id1`/`twin_id2` are temporarily twinless. Re-pair them with
                        // each other below; if one is already gone (dangling twin in a
                        // degenerate neighborhood) the other must still get a partner — never
                        // leave a surviving halfedge twinless (that would later make
                        // `remove_face` of its face abort mid-detach and break the chain).
                        let twin1_alive = self.halfedges.contains_key(twin_id1);
                        let twin2_alive = self.halfedges.contains_key(twin_id2);

                        // Sever any remaining one-sided back-references to the re-paired
                        // halfedges (see `weld_faces_at`): removing the two halves below
                        // touches no other twin pointer, so a pre-existing asymmetric
                        // pair can still reference `twin_id1`/`twin_id2` from elsewhere.
                        // Re-pair such a stray referencer with a fresh boundary half
                        // instead of leaving it twinless.
                        if twin1_alive
                            && let Some(t1) = self.halfedges.get(twin_id1).and_then(|h| h.twin)
                            && t1 != twin_id2
                            && let Some(t1_he) = self.halfedges.get(t1)
                            && t1_he.twin == Some(twin_id1)
                            && self
                                .pair_with_fresh_boundary_half(
                                    t1,
                                    self.halfedges[twin_id1].end_vertex,
                                )
                                .is_none()
                        {
                            error!("remove_neighboring_flaps: could not re-pair {t1:?}");
                        }
                        if twin2_alive
                            && let Some(t2) = self.halfedges.get(twin_id2).and_then(|h| h.twin)
                            && t2 != twin_id1
                            && let Some(t2_he) = self.halfedges.get(t2)
                            && t2_he.twin == Some(twin_id2)
                            && self
                                .pair_with_fresh_boundary_half(
                                    t2,
                                    self.halfedges[twin_id2].end_vertex,
                                )
                                .is_none()
                        {
                            error!("remove_neighboring_flaps: could not re-pair {t2:?}");
                        }

                        // The removed halves are opposite halves of one edge (`he1.end ==
                        // start_v_id2`), so the survivors take over as the vertex seeds:
                        // `twin_id1` is outgoing from `start_v_id2`, `twin_id2` from
                        // `start_v_id1`. When a survivor gets a fresh partner, that fresh
                        // half occupies the survivor's own start side instead.
                        match (twin1_alive, twin2_alive) {
                            (true, true) => {
                                self.halfedges.get_mut(twin_id1).unwrap().twin = Some(twin_id2);
                                self.halfedges.get_mut(twin_id2).unwrap().twin = Some(twin_id1);
                                #[cfg(feature = "instrumentation")]
                                crate::record_op_trace!(
                                    "flap removal: re-paired {twin_id1:?}<->{twin_id2:?}"
                                );
                                if let Some(v) = self.vertices.get_mut(start_v_id1) {
                                    v.outgoing_halfedge = Some(twin_id2);
                                }
                                if let Some(v) = self.vertices.get_mut(start_v_id2) {
                                    v.outgoing_halfedge = Some(twin_id1);
                                }
                            }
                            (true, false) => {
                                // `twin_id2` is gone; re-pair `twin_id1` with a fresh
                                // boundary half (outgoing from `start_v_id1`).
                                match self.pair_with_fresh_boundary_half(twin_id1, he1_end_v) {
                                    Some(fresh) => {
                                        if let Some(v) = self.vertices.get_mut(start_v_id1) {
                                            v.outgoing_halfedge = Some(fresh);
                                        }
                                        if let Some(v) = self.vertices.get_mut(start_v_id2) {
                                            v.outgoing_halfedge = Some(twin_id1);
                                        }
                                    }
                                    None => error!(
                                        "remove_neighboring_flaps: could not re-pair twin {twin_id1:?}"
                                    ),
                                }
                            }
                            (false, true) => {
                                match self.pair_with_fresh_boundary_half(twin_id2, he2_end_v) {
                                    Some(fresh) => {
                                        if let Some(v) = self.vertices.get_mut(start_v_id1) {
                                            v.outgoing_halfedge = Some(twin_id2);
                                        }
                                        if let Some(v) = self.vertices.get_mut(start_v_id2) {
                                            v.outgoing_halfedge = Some(fresh);
                                        }
                                    }
                                    None => error!(
                                        "remove_neighboring_flaps: could not re-pair twin {twin_id2:?}"
                                    ),
                                }
                            }
                            (false, false) => {
                                // Both twins of this edge are already gone, so nothing
                                // survives to take over as either endpoint's vertex seed.
                                // If a `outgoing_halfedge` seed still points at one of the
                                // just-removed halfedges, repair it so ring traversals
                                // don't start from a dead id.
                                self.reseed_outgoing_if_dead(start_v_id1);
                                self.reseed_outgoing_if_dead(start_v_id2);
                            }
                        }
                    }
                } else if halfedges_of_faces.len() == 1 {
                    tracing::debug!("Single orphaned halfedge after flap removal (boundary edge)");

                    if let Some(&single_he_id) = halfedges_of_faces.iter().next() {
                        if let Some(single_he) = self.halfedges.get(single_he_id) {
                            if single_he.twin.is_none() {
                                tracing::error!("Single orphaned halfedge found and has no twin");
                            }
                        } else {
                            tracing::error!("Single orphaned halfedge not found")
                        }
                    }
                }

                return Some(true);
            }
        }

        Some(false)
    }
}

#[cfg(test)]
mod tests {
    use crate::ops::AddOrGetEdge;

    use super::*;
    use glam::Vec3;

    #[test]
    fn test_remove_degenerate_faces() {
        crate::utils::get_tracing_subscriber();
        let mut meshgraph = MeshGraph::new();

        let center_v_id = meshgraph.add_vertex(Vec3::new(0.0, 0.0, 1.0));

        let v1_id = meshgraph.add_vertex(Vec3::new(-0.2, 0.0, 0.0));
        let he_c_1_id = meshgraph
            .add_or_get_edge(center_v_id, v1_id)
            .unwrap()
            .start_to_end_he_id;
        let v1p_id = meshgraph.add_vertex(Vec3::new(-0.2, 0.0, 0.0));
        let AddOrGetEdge {
            start_to_end_he_id: he_c_1p_id,
            twin_he_id: he_1p_c_id,
            ..
        } = meshgraph.add_or_get_edge(center_v_id, v1p_id).unwrap();

        let v2_id = meshgraph.add_vertex(Vec3::new(-1.0, 1.0, 0.0));
        let he_c_2_id = meshgraph
            .add_or_get_edge(center_v_id, v2_id)
            .unwrap()
            .start_to_end_he_id;

        let v3_id = meshgraph.add_vertex(Vec3::new(-1.0, -1.0, 0.0));
        let he_c_3_id = meshgraph
            .add_or_get_edge(center_v_id, v3_id)
            .unwrap()
            .start_to_end_he_id;

        let v4_id = meshgraph.add_vertex(Vec3::new(0.2, 0.0, 0.0));
        let he_c_4_id = meshgraph
            .add_or_get_edge(center_v_id, v4_id)
            .unwrap()
            .start_to_end_he_id;

        let v4p_id = meshgraph.add_vertex(Vec3::new(0.2, 0.0, 0.0));
        let AddOrGetEdge {
            start_to_end_he_id: he_c_4p_id,
            twin_he_id: he_4p_c_id,
            ..
        } = meshgraph.add_or_get_edge(center_v_id, v4p_id).unwrap();

        let v5_id = meshgraph.add_vertex(Vec3::new(1.0, -1.0, 0.0));
        let he_c_5_id = meshgraph
            .add_or_get_edge(center_v_id, v5_id)
            .unwrap()
            .start_to_end_he_id;

        let v6_id = meshgraph.add_vertex(Vec3::new(1.0, 1.0, 0.0));
        let he_c_6_id = meshgraph
            .add_or_get_edge(center_v_id, v6_id)
            .unwrap()
            .start_to_end_he_id;

        meshgraph
            .add_face_from_halfedges(he_c_1_id, he_c_2_id)
            .unwrap();
        meshgraph
            .add_face_from_halfedges(he_c_2_id, he_c_3_id)
            .unwrap();
        meshgraph
            .add_face_from_halfedges(he_c_3_id, he_c_1p_id)
            .unwrap();

        meshgraph
            .add_face_from_halfedges(he_c_1p_id, he_c_4p_id)
            .unwrap();

        meshgraph
            .add_face_from_halfedges(he_c_4p_id, he_c_5_id)
            .unwrap();
        meshgraph
            .add_face_from_halfedges(he_c_5_id, he_c_6_id)
            .unwrap();
        meshgraph
            .add_face_from_halfedges(he_c_6_id, he_c_4_id)
            .unwrap();

        meshgraph
            .add_face_from_halfedges(he_c_1_id, he_c_4_id)
            .unwrap();

        meshgraph.halfedges[he_c_1p_id].end_vertex = v1_id;
        meshgraph.halfedges[he_c_4p_id].end_vertex = v4_id;

        // already created from face above
        let AddOrGetEdge {
            start_to_end_he_id: he_1p_4p_id,
            twin_he_id: he_4p_1p_id,
            ..
        } = meshgraph.add_or_get_edge(v1p_id, v4p_id).unwrap();
        meshgraph.halfedges[he_1p_4p_id].end_vertex = v4_id;
        meshgraph.halfedges[he_4p_1p_id].end_vertex = v1_id;

        let AddOrGetEdge {
            start_to_end_he_id: he_3_1p_id,
            twin_he_id: he_1p_3_id,
            ..
        } = meshgraph.add_or_get_edge(v3_id, v1p_id).unwrap();
        meshgraph.halfedges[he_3_1p_id].end_vertex = v1_id;

        let AddOrGetEdge {
            start_to_end_he_id: he_4p_5_id,
            twin_he_id: he_5_4p_id,
            ..
        } = meshgraph.add_or_get_edge(v4p_id, v5_id).unwrap();
        meshgraph.halfedges[he_5_4p_id].end_vertex = v4_id;

        meshgraph.outgoing_halfedges[v1_id].push(he_1p_4p_id);
        meshgraph.outgoing_halfedges[v1_id].push(he_1p_c_id);
        meshgraph.outgoing_halfedges[v1_id].push(he_1p_3_id);
        meshgraph.outgoing_halfedges[v4_id].push(he_4p_1p_id);
        meshgraph.outgoing_halfedges[v4_id].push(he_4p_c_id);
        meshgraph.outgoing_halfedges[v4_id].push(he_4p_5_id);

        meshgraph.remove_only_vertex(v1p_id);
        meshgraph.remove_only_vertex(v4p_id);

        #[cfg(feature = "rerun")]
        meshgraph.log_rerun();

        let mut removed_vertices = vec![];
        let mut removed_halfedges = vec![];
        let mut removed_faces = vec![];

        meshgraph.remove_degenerate_faces(
            center_v_id,
            &mut removed_vertices,
            &mut removed_halfedges,
            &mut removed_faces,
        );

        #[cfg(feature = "rerun")]
        {
            meshgraph.log_rerun();
            crate::RR.flush_blocking().unwrap();
        }
    }

    #[test]
    fn test_remove_degenerate_edges() {
        crate::utils::get_tracing_subscriber();

        let mut meshgraph = MeshGraph::new();

        let center_v_id = meshgraph.add_vertex(Vec3::new(0.0, 0.0, 1.0));

        let v1_id = meshgraph.add_vertex(Vec3::new(0.0, 0.0, 0.0));
        let he_c_1_id = meshgraph
            .add_or_get_edge(center_v_id, v1_id)
            .unwrap()
            .start_to_end_he_id;

        let v1p_id = meshgraph.add_vertex(Vec3::new(0.0, 0.0, 0.0));
        let AddOrGetEdge {
            start_to_end_he_id: he_c_1p_id,
            twin_he_id: he_1p_c_id,
            ..
        } = meshgraph.add_or_get_edge(center_v_id, v1p_id).unwrap();

        let v2_id = meshgraph.add_vertex(Vec3::new(-1.0, 1.0, 0.0));
        let he_c_2_id = meshgraph
            .add_or_get_edge(center_v_id, v2_id)
            .unwrap()
            .start_to_end_he_id;

        let v3_id = meshgraph.add_vertex(Vec3::new(-1.0, -1.0, 0.0));
        let he_c_3_id = meshgraph
            .add_or_get_edge(center_v_id, v3_id)
            .unwrap()
            .start_to_end_he_id;

        let v5_id = meshgraph.add_vertex(Vec3::new(1.0, -1.0, 0.0));
        let he_c_5_id = meshgraph
            .add_or_get_edge(center_v_id, v5_id)
            .unwrap()
            .start_to_end_he_id;

        let v6_id = meshgraph.add_vertex(Vec3::new(1.0, 1.0, 0.0));
        let he_c_6_id = meshgraph
            .add_or_get_edge(center_v_id, v6_id)
            .unwrap()
            .start_to_end_he_id;

        meshgraph
            .add_face_from_halfedges(he_c_1_id, he_c_2_id)
            .unwrap();
        meshgraph
            .add_face_from_halfedges(he_c_2_id, he_c_3_id)
            .unwrap();
        meshgraph
            .add_face_from_halfedges(he_c_3_id, he_c_1p_id)
            .unwrap();

        meshgraph
            .add_face_from_halfedges(he_c_1p_id, he_c_5_id)
            .unwrap();

        meshgraph
            .add_face_from_halfedges(he_c_5_id, he_c_6_id)
            .unwrap();
        meshgraph
            .add_face_from_halfedges(he_c_6_id, he_c_1_id)
            .unwrap();

        meshgraph.halfedges[he_c_1p_id].end_vertex = v1_id;

        // already created from face above
        let AddOrGetEdge {
            start_to_end_he_id: he_3_1p_id,
            twin_he_id: he_1p_3_id,
            ..
        } = meshgraph.add_or_get_edge(v3_id, v1p_id).unwrap();
        meshgraph.halfedges[he_3_1p_id].end_vertex = v1_id;

        let AddOrGetEdge {
            start_to_end_he_id: he_1p_5_id,
            twin_he_id: he_5_1p_id,
            ..
        } = meshgraph.add_or_get_edge(v1p_id, v5_id).unwrap();
        meshgraph.halfedges[he_5_1p_id].end_vertex = v1_id;

        meshgraph.outgoing_halfedges[v1_id].push(he_1p_c_id);
        meshgraph.outgoing_halfedges[v1_id].push(he_1p_5_id);
        meshgraph.outgoing_halfedges[v1_id].push(he_1p_3_id);

        meshgraph.remove_only_vertex(v1p_id);

        #[cfg(feature = "rerun")]
        meshgraph.log_rerun();

        meshgraph.remove_degenerate_edges(center_v_id);

        #[cfg(feature = "rerun")]
        {
            meshgraph.log_rerun();
            crate::RR.flush_blocking().unwrap();
        }
    }

    #[test]
    fn test_remove_flap() {
        crate::utils::get_tracing_subscriber();

        let mut mesh_graph = MeshGraph::new();

        let v1 = mesh_graph.add_vertex(Vec3::new(0.0, 0.0, 0.0));
        let v2 = mesh_graph.add_vertex(Vec3::new(1.0, 0.0, 0.0));
        let v3 = mesh_graph.add_vertex(Vec3::new(0.0, 1.0, 0.0));

        let v4 = mesh_graph.add_vertex(Vec3::new(1.0, 1.0, 0.5));
        let v5 = mesh_graph.add_vertex(Vec3::new(1.0, 1.0, -0.5));

        let edge1 = mesh_graph.add_edge(v1, v2).unwrap();
        let edge2 = mesh_graph.add_edge(v2, v3).unwrap();
        let edge2_d = mesh_graph.add_edge(v2, v3).unwrap();

        mesh_graph
            .add_face_from_halfedges(edge1.start_to_end_he_id, edge2.start_to_end_he_id)
            .unwrap();
        mesh_graph
            .add_face_from_halfedges(edge2_d.twin_he_id, edge1.twin_he_id)
            .unwrap();

        mesh_graph
            .add_face_from_halfedge_and_vertex(edge2.twin_he_id, v4)
            .unwrap();
        mesh_graph
            .add_face_from_halfedge_and_vertex(edge2_d.start_to_end_he_id, v5)
            .unwrap();

        #[cfg(feature = "rerun")]
        mesh_graph.log_rerun();

        let mut removed_vertices = Vec::new();
        let mut removed_halfedges = Vec::new();
        let mut removed_faces = Vec::new();

        mesh_graph.remove_neighboring_flaps(
            v1,
            &mut removed_vertices,
            &mut removed_halfedges,
            &mut removed_faces,
        );

        #[cfg(feature = "rerun")]
        {
            mesh_graph.log_rerun();
            crate::RR.flush_blocking().unwrap();
        }

        assert_eq!(removed_vertices.len(), 1);
        assert_eq!(removed_halfedges.len(), 6);
        assert_eq!(removed_faces.len(), 2);
    }

    #[test]
    fn test_remove_double_flaps() {
        crate::utils::get_tracing_subscriber();

        let mut mesh_graph = MeshGraph::new();

        let v1 = mesh_graph.add_vertex(Vec3::new(0.0, 0.0, 0.0));
        let v2 = mesh_graph.add_vertex(Vec3::new(1.0, 0.0, 0.0));
        let v3 = mesh_graph.add_vertex(Vec3::new(0.0, 1.0, 0.0));
        let v4 = mesh_graph.add_vertex(Vec3::new(1.0, 1.0, 0.0));

        let v5 = mesh_graph.add_vertex(Vec3::new(0.0, -1.0, 0.0));

        let edge1 = mesh_graph.add_edge(v1, v2).unwrap();
        let edge1_d = mesh_graph.add_edge(v1, v2).unwrap();
        let edge2 = mesh_graph.add_edge(v2, v3).unwrap();
        let edge2_d = mesh_graph.add_edge(v2, v3).unwrap();

        mesh_graph
            .add_face_from_halfedges(edge1.start_to_end_he_id, edge2.start_to_end_he_id)
            .unwrap();
        mesh_graph
            .add_face_from_halfedge_and_vertex(edge2.twin_he_id, v4)
            .unwrap();

        mesh_graph
            .add_face_from_halfedge_and_vertex(edge2_d.start_to_end_he_id, v4)
            .unwrap();
        mesh_graph
            .add_face_from_halfedges(edge2_d.twin_he_id, edge1_d.twin_he_id)
            .unwrap();

        mesh_graph
            .add_face_from_halfedge_and_vertex(edge1_d.start_to_end_he_id, v5)
            .unwrap();

        #[cfg(feature = "rerun")]
        mesh_graph.log_rerun();

        let mut removed_vertices = Vec::new();
        let mut removed_halfedges = Vec::new();
        let mut removed_faces = Vec::new();

        mesh_graph.remove_neighboring_flaps(
            v1,
            &mut removed_vertices,
            &mut removed_halfedges,
            &mut removed_faces,
        );

        #[cfg(feature = "rerun")]
        {
            mesh_graph.log_rerun();
            crate::RR.flush_blocking().unwrap();
        }

        assert_eq!(removed_vertices.len(), 0);
        assert_eq!(removed_halfedges.len(), 6);
        assert_eq!(removed_faces.len(), 2);
    }
}