Fix #78916: unpredictable results with merge by distance

The previous implementation used KDTree duplicate search with
BLI_kdtree_3d_calc_duplicates_fast(). The survivor was always
one of the input vertices, not the centroid of the cluster.

This caused cases where merging a line of vertices did not
collapse to their average position, resulting in jagged loops.

Now vertices within the threshold are clustered, their centroid
is computed, and the chosen survivor is snapped to this centroid.
This ensures predictable and consistent merge results.

Ref !145851

source/blender/blenlib/BLI_kdtree_impl.h

@@ -68,6 +68,18 @@ int BLI_kdtree_nd_(calc_duplicates_fast)(const KDTree *tree,
                                          bool use_index_order,
                                          int *duplicates);
 
+/**
+ * Stable clustering and centroid computation to ensure consistent survivor selection.
+ *
+ * \param tree: A tree, all indices *must* be unique.
+ *
+ * \note ~1.1x-1.5x slower than `calc_duplicates_fast` depending on the distribution of points.
+ */
+int BLI_kdtree_nd_(calc_duplicates_and_center)(const KDTree *tree,
+                                               const float range,
+                                               int *duplicates,
+                                               float (*r_cluster_center)[KD_DIMS]);
+
 int BLI_kdtree_nd_(deduplicate)(KDTree *tree);
 
 /** Versions of find/range search that take a squared distance callback to support bias. */

source/blender/blenlib/intern/kdtree_impl.h

@@ -8,6 +8,8 @@
 
 #include "MEM_guardedalloc.h"
 
+#include "BLI_array.hh"
+#include "BLI_bit_vector.hh"
 #include "BLI_kdtree_impl.h"
 #include "BLI_math_base.h"
 #include "BLI_utildefines.h"
@@ -66,6 +68,20 @@ struct KDTree {
 /** \name Local Math API
  * \{ */
 
+static void add_vn_vn(float v0[KD_DIMS], const float v1[KD_DIMS])
+{
+  for (uint j = 0; j < KD_DIMS; j++) {
+    v0[j] += v1[j];
+  }
+}
+
+static void div_vn_fl(float v[KD_DIMS], const float f)
+{
+  for (uint j = 0; j < KD_DIMS; j++) {
+    v[j] /= f;
+  }
+}
+
 static void copy_vn_vn(float v0[KD_DIMS], const float v1[KD_DIMS])
 {
   for (uint j = 0; j < KD_DIMS; j++) {
@@ -934,6 +950,83 @@ int BLI_kdtree_nd_(calc_duplicates_fast)(const KDTree *tree,
 
 /** \} */
 
+/* -------------------------------------------------------------------- */
+/** \name BLI_kdtree_3d_calc_duplicates_and_center
+ * \{ */
+
+int BLI_kdtree_nd_(calc_duplicates_and_center)(const KDTree *tree,
+                                               const float range,
+                                               int *duplicates,
+                                               float (*r_cluster_center)[KD_DIMS])
+{
+  BLI_assert(tree->is_balanced);
+  if (UNLIKELY(tree->root == KD_NODE_UNSET)) {
+    return 0;
+  }
+
+  const uint nodes_len = tree->nodes_len;
+  blender::Array<int> index_to_node_index(tree->max_node_index + 1);
+  for (uint i = 0; i < nodes_len; i++) {
+    index_to_node_index[tree->nodes[i].index] = int(i);
+  }
+
+  blender::BitVector<> visited(tree->max_node_index + 1, false);
+  /* Could be inline, declare here to avoid re-allocation. */
+  blender::Vector<int> cluster;
+
+  int found = 0;
+
+  for (uint i = 0; i < nodes_len; i++) {
+    const int node_index = tree->nodes[i].index;
+    if ((duplicates[node_index] != -1) || visited[node_index]) {
+      continue;
+    }
+
+    const float *search_co = tree->nodes[index_to_node_index[node_index]].co;
+    BLI_assert(search_co != nullptr);
+    BLI_assert(cluster.is_empty());
+
+    visited[node_index].set();
+    auto accumulate_neighbors_fn = [&duplicates, &visited, &cluster](int neighbor_index,
+                                                                     const float * /*co*/,
+                                                                     float /*dist_sq*/) -> bool {
+      if ((duplicates[neighbor_index] == -1) && !visited[neighbor_index]) {
+        cluster.append(neighbor_index);
+        visited[neighbor_index].set();
+      }
+      return true;
+    };
+
+    BLI_kdtree_nd_(range_search_cb_cpp)(tree, search_co, range, accumulate_neighbors_fn);
+    if (cluster.is_empty()) {
+      continue;
+    }
+    cluster.append(node_index);
+
+    /* Compute centroid and choose survivor in one pass. */
+    float centroid[KD_DIMS] = {};
+    int survivor_index = node_index; /* Same as `cluster[0]`. */
+
+    for (int cluster_node_index : cluster) {
+      add_vn_vn(centroid, tree->nodes[index_to_node_index[cluster_node_index]].co);
+      survivor_index = std::min(cluster_node_index, survivor_index);
+    }
+
+    div_vn_fl(centroid, float(cluster.size()));
+    copy_vn_vn(r_cluster_center[survivor_index], centroid);
+
+    for (int cluster_node_index : cluster) {
+      duplicates[cluster_node_index] = survivor_index;
+    }
+    found += int(cluster.size()) - 1;
+    cluster.clear();
+  }
+
+  return found;
+}
+
+/** \} */
+
 /* -------------------------------------------------------------------- */
 /** \name BLI_kdtree_3d_deduplicate
  * \{ */

source/blender/bmesh/operators/bmo_removedoubles.cc

@@ -633,8 +633,13 @@ static int *bmesh_find_doubles_by_distance_impl(BMesh *bm,
                                                 const float dist,
                                                 const bool has_keep_vert)
 {
+  if (verts_len == 0) {
+    return nullptr;
+  }
+
   int *duplicates = MEM_malloc_arrayN<int>(verts_len, __func__);
   bool found_duplicates = false;
+  blender::Vector<blender::float3> survivor_cos(verts_len);
 
   KDTree_3d *tree = BLI_kdtree_3d_new(verts_len);
   for (int i = 0; i < verts_len; i++) {
@@ -648,13 +653,23 @@ static int *bmesh_find_doubles_by_distance_impl(BMesh *bm,
   }
 
   BLI_kdtree_3d_balance(tree);
-  found_duplicates = BLI_kdtree_3d_calc_duplicates_fast(tree, dist, false, duplicates) != 0;
+  found_duplicates = BLI_kdtree_3d_calc_duplicates_and_center(
+                         tree, dist, duplicates, (float(*)[3])survivor_cos.data()) != 0;
+
   BLI_kdtree_3d_free(tree);
 
   if (!found_duplicates) {
     MEM_freeN(duplicates);
     duplicates = nullptr;
   }
+  else {
+    for (int i = 0; i < verts_len; i++) {
+      if (duplicates[i] == i) {
+        copy_v3_v3(verts[i]->co, survivor_cos[i]);
+      }
+    }
+  }
+
   return duplicates;
 }