BMesh: New tool `BM_mesh_intersect_edges`

Along with the new utility `BM_vert_weld_linked_wire_edges_into_linked_faces`

1 files changed, 909 insertions, 0 deletions

diff --git a/source/blender/bmesh/tools/bmesh_intersect_edges.c b/source/blender/bmesh/tools/bmesh_intersect_edges.c
new file mode 100644
index 00000000000..04c6968bfa0
--- /dev/null
+++ b/source/blender/bmesh/tools/bmesh_intersect_edges.c
@@ -0,0 +1,909 @@
+/*
+ * This program is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU General Public License
+ * as published by the Free Software Foundation; either version 2
+ * of the License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program; if not, write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
+ *
+ * The Original Code is Copyright (C) 2019 Blender Foundation.
+ * All rights reserved.
+ */
+
+/** \file
+ * \ingroup bmesh
+ */
+
+#include "MEM_guardedalloc.h"
+
+#include "BLI_math.h"
+#include "BLI_sort.h"
+#include "BLI_stack.h"
+
+#include "BKE_bvhutils.h"
+
+#include "bmesh.h"
+
+#include "bmesh_intersect_edges.h" /* own include */
+
+#define KDOP_AXIS_LEN 14
+
+/* -------------------------------------------------------------------- */
+/** \name Weld Linked Wire Edges into Linked Faces
+ *
+ * Used with the merge vertices option.
+ * \{ */
+
+/* Callbacks for `BM_vert_pair_shared_face_cb` */
+
+struct EDBMSplitBestFaceData {
+  BMEdge **edgenet;
+  int edgenet_len;
+
+  /**
+   * Track the range of vertices in edgenet along the faces normal,
+   * find the lowest since it's most likely to be most co-planar with the face.
+   */
+  float best_face_range_on_normal_axis;
+  BMFace *r_best_face;
+};
+
+static bool bm_vert_pair_share_best_splittable_face_cb(BMFace *f,
+                                                       BMLoop *l_a,
+                                                       BMLoop *l_b,
+                                                       void *userdata)
+{
+  struct EDBMSplitBestFaceData *data = userdata;
+  float no[3];
+  copy_v3_v3(no, f->no);
+
+  float min = dot_v3v3(l_a->v->co, no);
+  float max = dot_v3v3(l_b->v->co, no);
+  if (min > max) {
+    SWAP(float, min, max);
+  }
+
+  BMVert *v_test = l_b->v;
+  BMEdge **e_iter = &data->edgenet[0];
+  int verts_len = data->edgenet_len - 1;
+  for (int i = verts_len; i--; e_iter++) {
+    v_test = BM_edge_other_vert(*e_iter, v_test);
+    if (!BM_face_point_inside_test(f, v_test->co)) {
+      return false;
+    }
+    float dot = dot_v3v3(v_test->co, no);
+    if (dot < min) {
+      min = dot;
+    }
+    if (dot > max) {
+      max = dot;
+    }
+  }
+
+  const float test_face_range_on_normal_axis = max - min;
+  if (test_face_range_on_normal_axis < data->best_face_range_on_normal_axis) {
+    data->best_face_range_on_normal_axis = test_face_range_on_normal_axis;
+    data->r_best_face = f;
+  }
+
+  return false;
+}
+
+/* find the best splittable face between the two vertices. */
+static bool bm_vert_pair_share_splittable_face_cb(BMFace *UNUSED(f),
+                                                  BMLoop *l_a,
+                                                  BMLoop *l_b,
+                                                  void *userdata)
+{
+  float(*data)[3] = userdata;
+  float *v_a_co = data[0];
+  float *v_a_b_dir = data[1];
+
+  float lambda;
+  if (isect_ray_seg_v3(v_a_co, v_a_b_dir, l_a->prev->v->co, l_a->next->v->co, &lambda)) {
+    if (IN_RANGE(lambda, 0.0f, 1.0f)) {
+      return true;
+    }
+    else if (isect_ray_seg_v3(v_a_co, v_a_b_dir, l_b->prev->v->co, l_b->next->v->co, &lambda)) {
+      return IN_RANGE(lambda, 0.0f, 1.0f);
+    }
+  }
+  return false;
+}
+
+void BM_vert_weld_linked_wire_edges_into_linked_faces(
+    BMesh *bm, BMVert *v, const float epsilon, BMEdge **r_edgenet[], int *r_edgenet_alloc_len)
+{
+  BMEdge **edgenet = *r_edgenet;
+  int edgenet_alloc_len = *r_edgenet_alloc_len;
+
+  BMIter iter;
+  BMEdge *e;
+  BM_ITER_ELEM (e, &iter, v, BM_EDGES_OF_VERT) {
+    int edgenet_len = 0;
+    BMVert *v_other = v;
+    while (BM_edge_is_wire(e)) {
+      if (edgenet_alloc_len == edgenet_len) {
+        edgenet_alloc_len = (edgenet_alloc_len + 1) * 2;
+        edgenet = MEM_reallocN(edgenet, (edgenet_alloc_len) * sizeof(*edgenet));
+      }
+      edgenet[edgenet_len++] = e;
+      v_other = BM_edge_other_vert(e, v_other);
+      if (v_other == v) {
+        /* Endless loop. */
+        break;
+      }
+
+      BMEdge *e_next = BM_DISK_EDGE_NEXT(e, v_other);
+      if (e_next == e) {
+        /* Vert is wire_endpoint. */
+        edgenet_len = 0;
+        break;
+      }
+
+      BMEdge *e_test = e_next;
+      while ((e_test = BM_DISK_EDGE_NEXT(e_test, v_other)) != e) {
+        if (e_test->l) {
+          /* Vert is linked to a face. */
+          goto l_break;
+        }
+      }
+
+      e = e_next;
+    }
+
+    BMLoop *dummy;
+    BMFace *best_face;
+
+  l_break:
+    if (edgenet_len == 0) {
+      /* Nothing to do. */
+      continue;
+    }
+    if (edgenet_len == 1) {
+      float data[2][3];
+      copy_v3_v3(data[0], v_other->co);
+      sub_v3_v3v3(data[1], v->co, data[0]);
+      best_face = BM_vert_pair_shared_face_cb(
+          v_other, v, true, bm_vert_pair_share_splittable_face_cb, &data, &dummy, &dummy);
+    }
+    else {
+      struct EDBMSplitBestFaceData data = {
+          .edgenet = edgenet,
+          .edgenet_len = edgenet_len,
+          .best_face_range_on_normal_axis = FLT_MAX,
+          .r_best_face = NULL,
+      };
+      BM_vert_pair_shared_face_cb(
+          v_other, v, true, bm_vert_pair_share_best_splittable_face_cb, &data, &dummy, &dummy);
+
+      if (data.r_best_face) {
+        float no[3], min = FLT_MAX, max = -FLT_MAX;
+        copy_v3_v3(no, data.r_best_face->no);
+        BMVert *v_test;
+        BMIter f_iter;
+        BM_ITER_ELEM (v_test, &f_iter, data.r_best_face, BM_VERTS_OF_FACE) {
+          float dot = dot_v3v3(v_test->co, no);
+          if (dot < min) {
+            min = dot;
+          }
+          if (dot > max) {
+            max = dot;
+          }
+        }
+        float range = max - min + 2 * epsilon;
+        if (range < data.best_face_range_on_normal_axis) {
+          data.r_best_face = NULL;
+        }
+      }
+      best_face = data.r_best_face;
+    }
+
+    if (best_face) {
+      BM_face_split_edgenet(bm, best_face, edgenet, edgenet_len, NULL, NULL);
+    }
+  }
+
+  *r_edgenet = edgenet;
+  *r_edgenet_alloc_len = edgenet_alloc_len;
+}
+
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Auto-Merge & Split Selection
+ *
+ * Used after transform operations.
+ * \{ */
+
+struct EDBMSplitElem {
+  union {
+    BMElem *elem;
+    BMVert *vert;
+    struct {
+      BMEdge *edge;
+      float lambda;
+    };
+  };
+};
+
+/* -------------------------------------------------------------------- */
+/* Overlap Callbacks */
+
+struct EDBMSplitData {
+  BMesh *bm;
+  BLI_Stack *pair_stack;
+  int cut_edges_a_len;
+  int cut_edges_b_len;
+  float dist_sq;
+  float dist_sq_sq;
+};
+
+/* Utils */
+
+static void bm_vert_pair_elem_setup_ex(BMVert *v,
+                                       float edge_index,
+                                       struct EDBMSplitElem *r_pair_elem)
+{
+  BLI_assert(v->head.index == -1);
+  v->head.index = edge_index;
+  r_pair_elem->vert = v;
+}
+
+static void bm_edge_pair_elem_setup(BMEdge *e,
+                                    float lambda,
+                                    int *r_data_cut_edges_len,
+                                    struct EDBMSplitElem *r_pair_elem)
+{
+  r_pair_elem->edge = e;
+  r_pair_elem->lambda = lambda;
+
+  e->head.index++;
+  /* Obs: Check Multithread. */
+  if (BM_elem_flag_test(e, BM_ELEM_TAG)) {
+    BM_elem_flag_disable(e, BM_ELEM_TAG);
+    (*r_data_cut_edges_len)++;
+  }
+}
+
+/* Util for Vert x Edge and Edge x Edge callbacks */
+static bool bm_vertxedge_isect_impl_ex(BMVert *v,
+                                       BMEdge *e,
+                                       int edge_index,
+                                       float co[3],
+                                       float dir[3],
+                                       float lambda,
+                                       float data_dist_sq,
+                                       int *data_cut_edges_len,
+                                       struct EDBMSplitElem r_pair[2])
+{
+  BLI_assert(v->head.index == -1);
+
+  BMVert *e_v;
+  float dist_sq_vert_factor;
+
+  if (lambda < 0.5f) {
+    e_v = e->v1;
+    dist_sq_vert_factor = lambda;
+  }
+  else {
+    e_v = e->v2;
+    dist_sq_vert_factor = 1.0f - lambda;
+  }
+
+  if (v != e_v) {
+    CLAMP(lambda, 0.0f, 1.0f);
+
+    float near[3];
+    madd_v3_v3v3fl(near, co, dir, lambda);
+
+    float dist_sq = len_squared_v3v3(v->co, near);
+    if (dist_sq < data_dist_sq) {
+      float dist_sq_vert = SQUARE(dist_sq_vert_factor) * len_squared_v3(dir);
+      if (dist_sq_vert < data_dist_sq) {
+        if (e_v->head.index != -1) {
+          /* Vertex already has an intersection. */
+          return false;
+        }
+
+        bm_vert_pair_elem_setup_ex(e_v, -2, &r_pair[1]);
+      }
+      else {
+        bm_edge_pair_elem_setup(e, lambda, data_cut_edges_len, &r_pair[1]);
+      }
+
+      bm_vert_pair_elem_setup_ex(v, edge_index, &r_pair[0]);
+      return true;
+    }
+  }
+
+  return false;
+}
+
+/* Vertex x Vertex Callback */
+
+static bool bm_vertxvert_isect_cb(void *userdata, int index_a, int index_b, int UNUSED(thread))
+{
+  struct EDBMSplitData *data = userdata;
+  BMVert *v_a = BM_vert_at_index(data->bm, index_a);
+  BMVert *v_b = BM_vert_at_index(data->bm, index_b);
+
+  struct EDBMSplitElem *pair = BLI_stack_push_r(data->pair_stack);
+
+  BLI_assert(v_a->head.index == -1);
+
+  /* Set index -2 for sure that it will not repeat keys in `targetmap`. */
+  bm_vert_pair_elem_setup_ex(v_a, -2, &pair[0]);
+  bm_vert_pair_elem_setup_ex(v_b, -1, &pair[1]);
+
+  return true;
+}
+
+/* Vertex x Edge and Edge x Vertex Callbacks */
+
+static int bm_vertxedge_isect_impl(BMesh *bm,
+                                   int vert_index,
+                                   int edge_index,
+                                   float data_dist_sq,
+                                   int *data_cut_edges_len,
+                                   struct EDBMSplitElem r_pair[2])
+{
+  BMVert *v = BM_vert_at_index(bm, vert_index);
+  BMEdge *e = BM_edge_at_index(bm, edge_index);
+
+  if (v->head.index != -1) {
+    /* Only one vertex per edge. */
+    return false;
+  }
+
+  float co[3], dir[3], lambda;
+  copy_v3_v3(co, e->v1->co);
+  sub_v3_v3v3(dir, e->v2->co, co);
+  lambda = ray_point_factor_v3_ex(v->co, co, dir, 0.0f, -1.0f);
+
+  return bm_vertxedge_isect_impl_ex(
+      v, e, edge_index, co, dir, lambda, data_dist_sq, data_cut_edges_len, r_pair);
+}
+
+static bool bm_vertxedge_isect_cb(void *userdata, int index_a, int index_b, int UNUSED(thread))
+{
+  struct EDBMSplitData *data = userdata;
+  struct EDBMSplitElem pair_tmp[2];
+  if (bm_vertxedge_isect_impl(
+          data->bm, index_a, index_b, data->dist_sq, &data->cut_edges_b_len, pair_tmp)) {
+    struct EDBMSplitElem *pair = BLI_stack_push_r(data->pair_stack);
+    pair[0] = pair_tmp[0];
+    pair[1] = pair_tmp[1];
+
+    return true;
+  }
+
+  return false;
+}
+
+static bool bm_edgexvert_isect_cb(void *userdata, int index_a, int index_b, int UNUSED(thread))
+{
+  struct EDBMSplitData *data = userdata;
+  struct EDBMSplitElem pair_tmp[2];
+  if (bm_vertxedge_isect_impl(
+          data->bm, index_b, index_a, data->dist_sq, &data->cut_edges_a_len, pair_tmp)) {
+    struct EDBMSplitElem *pair = BLI_stack_push_r(data->pair_stack);
+    pair[0] = pair_tmp[1];
+    pair[1] = pair_tmp[0];
+
+    return true;
+  }
+
+  return false;
+}
+
+/* Edge x Edge Callbacks */
+
+static void bm_edgexedge_isect_impl(struct EDBMSplitData *data,
+                                    int index_a,
+                                    int index_b,
+                                    BMEdge *e_a,
+                                    BMEdge *e_b,
+                                    float co_a[3],
+                                    float dir_a[3],
+                                    float co_b[3],
+                                    float dir_b[3],
+                                    float lambda_a,
+                                    float lambda_b)
+{
+  float dist_sq_va_factor, dist_sq_vb_factor;
+  BMVert *e_a_v, *e_b_v;
+  if (lambda_a < 0.5f) {
+    e_a_v = e_a->v1;
+    dist_sq_va_factor = lambda_a;
+  }
+  else {
+    e_a_v = e_a->v2;
+    dist_sq_va_factor = 1.0f - lambda_a;
+  }
+
+  if (lambda_b < 0.5f) {
+    e_b_v = e_b->v1;
+    dist_sq_vb_factor = lambda_b;
+  }
+  else {
+    e_b_v = e_b->v2;
+    dist_sq_vb_factor = 1.0f - lambda_b;
+  }
+
+  if (e_a_v != e_b_v) {
+    CLAMP(lambda_a, 0.0f, 1.0f);
+    CLAMP(lambda_b, 0.0f, 1.0f);
+
+    float near_a[3], near_b[3];
+    madd_v3_v3v3fl(near_a, co_a, dir_a, lambda_a);
+    madd_v3_v3v3fl(near_b, co_b, dir_b, lambda_b);
+
+    float dist_sq = len_squared_v3v3(near_a, near_b);
+    if (dist_sq < data->dist_sq) {
+      struct EDBMSplitElem pair_tmp[2];
+
+      float dist_sq_va = SQUARE(dist_sq_va_factor) * len_squared_v3(dir_a);
+      float dist_sq_vb = SQUARE(dist_sq_vb_factor) * len_squared_v3(dir_b);
+
+      if (dist_sq_va < data->dist_sq) {
+        if (e_a_v->head.index != -1) {
+          /* Only one vertex per edge. */
+          return;
+        }
+        bm_vert_pair_elem_setup_ex(e_a_v, index_b, &pair_tmp[0]);
+      }
+
+      if (dist_sq_vb < data->dist_sq) {
+        if (e_b_v->head.index != -1) {
+          /* Only one vertex per edge. */
+          return;
+        }
+        bm_vert_pair_elem_setup_ex(e_b_v, index_a, &pair_tmp[1]);
+      }
+      else {
+        bm_edge_pair_elem_setup(e_b, lambda_b, &data->cut_edges_b_len, &pair_tmp[1]);
+      }
+
+      /* Don't setup edges before a return. */
+      if (dist_sq_va >= data->dist_sq) {
+        bm_edge_pair_elem_setup(e_a, lambda_a, &data->cut_edges_a_len, &pair_tmp[0]);
+      }
+
+      struct EDBMSplitElem *pair = BLI_stack_push_r(data->pair_stack);
+      pair[0] = pair_tmp[0];
+      pair[1] = pair_tmp[1];
+    }
+  }
+}
+
+static bool bm_edgexedge_isect_cb(void *userdata, int index_a, int index_b, int UNUSED(thread))
+{
+  bool ret = false;
+  struct EDBMSplitData *data = userdata;
+  BMEdge *e_a = BM_edge_at_index(data->bm, index_a);
+  BMEdge *e_b = BM_edge_at_index(data->bm, index_b);
+
+  float co_a[3], dir_a[3], co_b[3], dir_b[3];
+  copy_v3_v3(co_a, e_a->v1->co);
+  sub_v3_v3v3(dir_a, e_a->v2->co, co_a);
+
+  copy_v3_v3(co_b, e_b->v1->co);
+  sub_v3_v3v3(dir_b, e_b->v2->co, co_b);
+
+  float lambda_a, lambda_b;
+  /* Using with dist^4 as `epsilon` is not the best solution, but it fits in most cases. */
+  if (isect_ray_ray_epsilon_v3(co_a, dir_a, co_b, dir_b, data->dist_sq_sq, &lambda_a, &lambda_b)) {
+    if (ELEM(index_b, e_a->v1->head.index, e_a->v2->head.index) ||
+        ELEM(index_a, e_b->v1->head.index, e_b->v2->head.index)) {
+      return ret;
+    }
+
+    /* Edge x Edge returns always false. */
+    bm_edgexedge_isect_impl(
+        data, index_a, index_b, e_a, e_b, co_a, dir_a, co_b, dir_b, lambda_a, lambda_b);
+  }
+  else {
+    /* Parallel */
+    struct EDBMSplitElem pair_tmp[2];
+    float vec[3], len_sq_a, len_sq_b, lambda;
+    sub_v3_v3v3(vec, co_b, co_a);
+    len_sq_a = len_squared_v3(dir_a);
+    len_sq_b = len_squared_v3(dir_b);
+
+    if (!ELEM(e_b->v1, e_a->v1, e_a->v2) && e_b->v1->head.index == -1) {
+      lambda = dot_v3v3(vec, dir_a) / len_sq_a;
+      if (bm_vertxedge_isect_impl_ex(e_b->v1,
+                                     e_a,
+                                     index_a,
+                                     co_a,
+                                     dir_a,
+                                     lambda,
+                                     data->dist_sq,
+                                     &data->cut_edges_a_len,
+                                     pair_tmp)) {
+        struct EDBMSplitElem *pair = BLI_stack_push_r(data->pair_stack);
+        pair[0] = pair_tmp[1];
+        pair[1] = pair_tmp[0];
+        ret |= true;
+      }
+    }
+
+    if (!ELEM(e_a->v1, e_b->v1, e_b->v2) && e_a->v1->head.index == -1) {
+      lambda = -dot_v3v3(vec, dir_b) / len_sq_b;
+      if (bm_vertxedge_isect_impl_ex(e_a->v1,
+                                     e_b,
+                                     index_b,
+                                     co_b,
+                                     dir_b,
+                                     lambda,
+                                     data->dist_sq,
+                                     &data->cut_edges_b_len,
+                                     pair_tmp)) {
+        struct EDBMSplitElem *pair = BLI_stack_push_r(data->pair_stack);
+        pair[0] = pair_tmp[0];
+        pair[1] = pair_tmp[1];
+        ret |= true;
+      }
+    }
+
+    add_v3_v3(vec, dir_b);
+    if (!ELEM(e_b->v2, e_a->v1, e_a->v2) && e_b->v2->head.index == -1) {
+      lambda = dot_v3v3(vec, dir_a) / len_sq_a;
+      if (bm_vertxedge_isect_impl_ex(e_b->v2,
+                                     e_a,
+                                     index_a,
+                                     co_a,
+                                     dir_a,
+                                     lambda,
+                                     data->dist_sq,
+                                     &data->cut_edges_a_len,
+                                     pair_tmp)) {
+        struct EDBMSplitElem *pair = BLI_stack_push_r(data->pair_stack);
+        pair[0] = pair_tmp[1];
+        pair[1] = pair_tmp[0];
+        ret |= true;
+      }
+    }
+
+    sub_v3_v3(vec, dir_a);
+    if (!ELEM(e_a->v2, e_b->v1, e_b->v2) && e_a->v2->head.index == -1) {
+      lambda = 1.0f - dot_v3v3(vec, dir_b) / len_sq_b;
+      if (bm_vertxedge_isect_impl_ex(e_a->v2,
+                                     e_b,
+                                     index_b,
+                                     co_b,
+                                     dir_b,
+                                     lambda,
+                                     data->dist_sq,
+                                     &data->cut_edges_b_len,
+                                     pair_tmp)) {
+        struct EDBMSplitElem *pair = BLI_stack_push_r(data->pair_stack);
+        pair[0] = pair_tmp[0];
+        pair[1] = pair_tmp[1];
+        ret |= true;
+      }
+    }
+  }
+
+  return ret;
+}
+
+/* -------------------------------------------------------------------- */
+/* BVHTree Overlap Function */
+
+static void bvhtree_overlap_thread_safe(const BVHTree *tree1,
+                                        const BVHTree *tree2,
+                                        BVHTree_OverlapCallback callback,
+                                        void *userdata)
+{
+  BLI_bvhtree_overlap_ex(tree1, tree2, NULL, callback, userdata, BVH_OVERLAP_BREAK_ON_FIRST);
+}
+
+/* -------------------------------------------------------------------- */
+/* Callbacks for `BLI_qsort_r` */
+
+static int sort_cmp_by_lambda_a_cb(const void *index1_v, const void *index2_v, void *keys_v)
+{
+  const struct EDBMSplitElem(*pair_array)[2] = keys_v;
+  const int index1 = *(int *)index1_v;
+  const int index2 = *(int *)index2_v;
+
+  if (pair_array[index1][0].lambda > pair_array[index2][0].lambda) {
+    return 1;
+  }
+  else {
+    return -1;
+  }
+}
+
+static int sort_cmp_by_lambda_b_cb(const void *index1_v, const void *index2_v, void *keys_v)
+{
+  const struct EDBMSplitElem(*pair_array)[2] = keys_v;
+  const int index1 = *(int *)index1_v;
+  const int index2 = *(int *)index2_v;
+
+  if (pair_array[index1][1].lambda > pair_array[index2][1].lambda) {
+    return 1;
+  }
+  else {
+    return -1;
+  }
+}
+
+/* -------------------------------------------------------------------- */
+/* Main API */
+
+bool BM_mesh_intersect_edges(BMesh *bm, const char hflag, const float dist, GHash *r_targetmap)
+{
+  bool ok = false;
+
+  BMIter iter;
+  BMVert *v;
+  BMEdge *e;
+  int i;
+
+  /* Store all intersections in this array. */
+  struct EDBMSplitElem(*pair_iter)[2], (*pair_array)[2] = NULL;
+  BLI_Stack *pair_stack = BLI_stack_new(sizeof(*pair_array), __func__);
+  int pair_len = 0;
+
+  float dist_sq = SQUARE(dist);
+  struct EDBMSplitData data = {
+      .bm = bm,
+      .pair_stack = pair_stack,
+      .cut_edges_a_len = 0,
+      .cut_edges_b_len = 0,
+      .dist_sq = dist_sq,
+      .dist_sq_sq = SQUARE(dist_sq),
+  };
+
+  /* tag and count the verts to be tested. */
+  int verts_act_len = 0, verts_remain_len = 0;
+  int loose_verts_act_len = 0, loose_verts_remain_len = 0;
+  BM_ITER_MESH (v, &iter, bm, BM_VERTS_OF_MESH) {
+    if (BM_elem_flag_test(v, hflag)) {
+      BM_elem_flag_enable(v, BM_ELEM_TAG);
+      v->head.index = -1;
+      verts_act_len++;
+      if (!v->e) {
+        loose_verts_act_len++;
+      }
+    }
+    else {
+      BM_elem_flag_disable(v, BM_ELEM_TAG);
+      if (!BM_elem_flag_test(v, BM_ELEM_HIDDEN)) {
+        v->head.index = -1;
+        verts_remain_len++;
+        if (!v->e) {
+          loose_verts_remain_len++;
+        }
+      }
+    }
+  }
+  bm->elem_index_dirty |= BM_VERT;
+
+  /* Start the creation of BVHTrees. */
+  BVHTree *tree_loose_verts_act = NULL, *tree_loose_verts_remain = NULL;
+  if (loose_verts_act_len) {
+    tree_loose_verts_act = BLI_bvhtree_new(loose_verts_act_len, dist, 2, KDOP_AXIS_LEN);
+  }
+  if (loose_verts_remain_len) {
+    tree_loose_verts_remain = BLI_bvhtree_new(loose_verts_remain_len, 0.0f, 2, KDOP_AXIS_LEN);
+  }
+
+  if (tree_loose_verts_act || tree_loose_verts_remain) {
+    BM_ITER_MESH_INDEX (v, &iter, bm, BM_VERTS_OF_MESH, i) {
+      if (BM_elem_flag_test(v, BM_ELEM_TAG)) {
+        if (tree_loose_verts_act && !v->e) {
+          BLI_bvhtree_insert(tree_loose_verts_act, i, v->co, 1);
+        }
+      }
+      else if (tree_loose_verts_remain && !v->e && !BM_elem_flag_test(v, BM_ELEM_HIDDEN)) {
+        BLI_bvhtree_insert(tree_loose_verts_remain, i, v->co, 1);
+      }
+    }
+    if (tree_loose_verts_act) {
+      BLI_bvhtree_balance(tree_loose_verts_act);
+    }
+
+    if (tree_loose_verts_remain) {
+      BLI_bvhtree_balance(tree_loose_verts_remain);
+    }
+
+    if (tree_loose_verts_act && tree_loose_verts_remain) {
+      /* First pair search. */
+      bvhtree_overlap_thread_safe(
+          tree_loose_verts_act, tree_loose_verts_remain, bm_vertxvert_isect_cb, &data);
+    }
+  }
+
+  /* Tag and count the edges. */
+  int edges_act_len = 0, edges_remain_len = 0;
+  BM_ITER_MESH (e, &iter, bm, BM_EDGES_OF_MESH) {
+    if (BM_elem_flag_test(e->v1, BM_ELEM_TAG) || BM_elem_flag_test(e->v2, BM_ELEM_TAG)) {
+      BM_elem_flag_enable(e, BM_ELEM_TAG);
+      edges_act_len++;
+    }
+    else {
+      BM_elem_flag_disable(e, BM_ELEM_TAG);
+      if (!BM_elem_flag_test(e, BM_ELEM_HIDDEN)) {
+        edges_remain_len++;
+      }
+    }
+  }
+
+  if (edges_remain_len) {
+    BVHTree *tree_edges_act = NULL, *tree_edges_remain = NULL;
+    tree_edges_remain = BLI_bvhtree_new(edges_remain_len, 0.0f, 2, KDOP_AXIS_LEN);
+    if (edges_act_len) {
+      tree_edges_act = BLI_bvhtree_new(edges_act_len, dist, 2, KDOP_AXIS_LEN);
+    }
+
+    BM_ITER_MESH_INDEX (e, &iter, bm, BM_EDGES_OF_MESH, i) {
+      float co[2][3];
+      if (BM_elem_flag_test(e, BM_ELEM_TAG)) {
+        if (tree_edges_act) {
+          e->head.index = 0;
+          copy_v3_v3(co[0], e->v1->co);
+          copy_v3_v3(co[1], e->v2->co);
+          BLI_bvhtree_insert(tree_edges_act, i, co[0], 2);
+        }
+      }
+      else if (!BM_elem_flag_test(e, BM_ELEM_HIDDEN)) {
+        /* Tag used in the overlap callbacks. */
+        BM_elem_flag_enable(e, BM_ELEM_TAG);
+        e->head.index = 0;
+        copy_v3_v3(co[0], e->v1->co);
+        copy_v3_v3(co[1], e->v2->co);
+        BLI_bvhtree_insert(tree_edges_remain, i, co[0], 2);
+      }
+    }
+    /* Use `e->head.index` to count intersections. */
+    bm->elem_index_dirty |= BM_EDGE;
+
+    BLI_bvhtree_balance(tree_edges_remain);
+    if (tree_edges_act) {
+      BLI_bvhtree_balance(tree_edges_act);
+    }
+
+    if (tree_edges_act) {
+      /* Edge x Edge */
+      bvhtree_overlap_thread_safe(tree_edges_act, tree_edges_remain, bm_edgexedge_isect_cb, &data);
+
+      if (tree_loose_verts_remain) {
+        /* Edge x Vert */
+        bvhtree_overlap_thread_safe(
+            tree_edges_act, tree_loose_verts_remain, bm_edgexvert_isect_cb, &data);
+      }
+
+      BLI_bvhtree_free(tree_edges_act);
+    }
+
+    if (tree_loose_verts_act) {
+      /* Vert x Edge */
+      bvhtree_overlap_thread_safe(
+          tree_loose_verts_act, tree_edges_remain, bm_vertxedge_isect_cb, &data);
+    }
+
+    BLI_bvhtree_free(tree_edges_remain);
+
+    pair_len = BLI_stack_count(pair_stack);
+    if (pair_len) {
+      pair_array = MEM_mallocN(sizeof(*pair_array) * pair_len, __func__);
+      BLI_stack_pop_n_reverse(pair_stack, pair_array, pair_len);
+
+      /* Map intersections per edge. */
+      union {
+        struct {
+          int cuts_len;
+          int cuts_index[];
+        };
+        int as_int[0];
+      } * e_map_iter, *e_map;
+
+      size_t e_map_size = (max_ii(data.cut_edges_a_len, data.cut_edges_b_len) * sizeof(*e_map)) +
+                          (pair_len * sizeof(*(e_map->cuts_index)));
+
+      e_map = MEM_mallocN(e_map_size, __func__);
+
+      /* Convert every pair to Vert x Vert. */
+      for (int pair = 0; pair < 2; pair++) {
+        int map_len = 0;
+        pair_iter = &pair_array[0];
+        for (i = 0; i < pair_len; i++, pair_iter++) {
+          if ((*pair_iter)[pair].elem->head.htype != BM_EDGE) {
+            /* Take the opportunity to set all vert indices to -1 again. */
+            (*pair_iter)[pair].elem->head.index = -1;
+            continue;
+          }
+          e = (*pair_iter)[pair].edge;
+          if (!BM_elem_flag_test(e, BM_ELEM_TAG)) {
+            BM_elem_flag_enable(e, BM_ELEM_TAG);
+            int e_cuts_len = e->head.index;
+
+            e_map_iter = (void *)&e_map->as_int[map_len];
+            e_map_iter->cuts_len = e_cuts_len;
+            e_map_iter->cuts_index[0] = i;
+
+            /* Use `e->head.index` to indicate which slot to fill with the `cut` index. */
+            e->head.index = map_len + 1;
+            map_len += 1 + e_cuts_len;
+          }
+          else {
+            e_map->as_int[++e->head.index] = i;
+          }
+        }
+
+        /* Split Edges A to set all Vert x Edge. */
+        for (i = 0; i < map_len;
+             e_map_iter = (void *)&e_map->as_int[i], i += 1 + e_map_iter->cuts_len) {
+
+          /* sort by lambda. */
+          BLI_qsort_r(e_map_iter->cuts_index,
+                      e_map_iter->cuts_len,
+                      sizeof(*(e_map->cuts_index)),
+                      pair == 0 ? sort_cmp_by_lambda_a_cb : sort_cmp_by_lambda_b_cb,
+                      pair_array);
+
+          float lambda, lambda_prev = 0.0f;
+          for (int j = 0; j < e_map_iter->cuts_len; j++) {
+            struct EDBMSplitElem *pair_elem = &pair_array[e_map_iter->cuts_index[j]][pair];
+            lambda = (pair_elem->lambda - lambda_prev) / (1.0f - lambda_prev);
+            lambda_prev = pair_elem->lambda;
+            e = pair_elem->edge;
+
+            BMVert *v_new = BM_edge_split(bm, e, e->v1, NULL, lambda);
+            v_new->head.index = -1;
+            pair_elem->vert = v_new;
+          }
+        }
+      }
+
+      MEM_freeN(e_map);
+    }
+  }
+
+  BLI_bvhtree_free(tree_loose_verts_act);
+  BLI_bvhtree_free(tree_loose_verts_remain);
+
+  if (r_targetmap) {
+    if (pair_array == NULL) {
+      pair_len = BLI_stack_count(pair_stack);
+      if (pair_len) {
+        pair_array = MEM_mallocN(sizeof(*pair_array) * pair_len, __func__);
+        BLI_stack_pop_n_reverse(pair_stack, pair_array, pair_len);
+      }
+    }
+
+    if (pair_array) {
+      /* Organize the vertices in the order they will be merged. */
+      pair_iter = &pair_array[0];
+      for (i = 0; i < pair_len; i++, pair_iter++) {
+        BLI_assert((*pair_iter)[0].elem->head.htype == BM_VERT);
+        BLI_assert((*pair_iter)[1].elem->head.htype == BM_VERT);
+        BLI_assert((*pair_iter)[0].elem != (*pair_iter)[1].elem);
+
+        BLI_ghash_insert(r_targetmap, (*pair_iter)[0].vert, (*pair_iter)[1].vert);
+      }
+
+      ok = true;
+    }
+  }
+
+  BLI_stack_free(pair_stack);
+  if (pair_array) {
+    MEM_freeN(pair_array);
+  }
+
+  return ok;
+}
+
+/** \} */