diff --git a/source/blender/blenkernel/BKE_paint.h b/source/blender/blenkernel/BKE_paint.h
index 0a35658464a..a8b58993864 100644
--- a/source/blender/blenkernel/BKE_paint.h
+++ b/source/blender/blenkernel/BKE_paint.h
@@ -473,6 +473,12 @@ typedef struct SculptSession {
   struct MeshElemMap *pmap;
   int *pmap_mem;
 
+  struct MeshElemMap *epmap;
+  int *epmap_mem;
+
+  struct MeshElemMap *vemap;
+  int *vemap_mem;
+
   /* Mesh Face Sets */
   /* Total number of polys of the base mesh. */
   int totfaces;
diff --git a/source/blender/blenkernel/intern/paint.c b/source/blender/blenkernel/intern/paint.c
index 4eecf3a3a87..5e4535ded56 100644
--- a/source/blender/blenkernel/intern/paint.c
+++ b/source/blender/blenkernel/intern/paint.c
@@ -1416,6 +1416,12 @@ static void sculptsession_free_pbvh(Object *object)
   MEM_SAFE_FREE(ss->pmap);
   MEM_SAFE_FREE(ss->pmap_mem);
 
+  MEM_SAFE_FREE(ss->epmap);
+  MEM_SAFE_FREE(ss->epmap_mem);
+
+  MEM_SAFE_FREE(ss->vemap);
+  MEM_SAFE_FREE(ss->vemap_mem);
+
   MEM_SAFE_FREE(ss->persistent_base);
 
   MEM_SAFE_FREE(ss->preview_vert_index_list);
@@ -1465,6 +1471,10 @@ void BKE_sculptsession_free(Object *ob)
 
     MEM_SAFE_FREE(ss->pmap);
     MEM_SAFE_FREE(ss->pmap_mem);
+    MEM_SAFE_FREE(ss->epmap);
+    MEM_SAFE_FREE(ss->epmap_mem);
+    MEM_SAFE_FREE(ss->vemap);
+    MEM_SAFE_FREE(ss->vemap_mem);
     if (ss->bm_log) {
       BM_log_free(ss->bm_log);
     }
diff --git a/source/blender/editors/sculpt_paint/sculpt.c b/source/blender/editors/sculpt_paint/sculpt.c
index 5070dba14ea..d46c23a9aa5 100644
--- a/source/blender/editors/sculpt_paint/sculpt.c
+++ b/source/blender/editors/sculpt_paint/sculpt.c
@@ -29,6 +29,7 @@
 #include "BLI_ghash.h"
 #include "BLI_gsqueue.h"
 #include "BLI_hash.h"
+#include "BLI_linklist_stack.h"
 #include "BLI_math.h"
 #include "BLI_math_color_blend.h"
 #include "BLI_task.h"
@@ -3961,33 +3962,36 @@ static void do_grab_brush_task_cb_ex(void *__restrict userdata,
   {
     SCULPT_orig_vert_data_update(&orig_data, &vd);
 
-    if (sculpt_brush_test_sq_fn(&test, orig_data.co)) {
-      float fade = bstrength * SCULPT_brush_strength_factor(ss,
-                                                            brush,
-                                                            orig_data.co,
-                                                            sqrtf(test.dist),
-                                                            orig_data.no,
-                                                            NULL,
-                                                            vd.mask ? *vd.mask : 0.0f,
-                                                            vd.index,
-                                                            thread_id);
-
-      if (grab_silhouette) {
-        float silhouette_test_dir[3];
-        normalize_v3_v3(silhouette_test_dir, grab_delta);
-        if (dot_v3v3(ss->cache->initial_normal, ss->cache->grab_delta_symmetry) < 0.0f) {
-          mul_v3_fl(silhouette_test_dir, -1.0f);
-        }
-        float vno[3];
-        normal_short_to_float_v3(vno, orig_data.no);
-        fade *= max_ff(dot_v3v3(vno, silhouette_test_dir), 0.0f);
-      }
+    if (!sculpt_brush_test_sq_fn(&test, orig_data.co)) {
+      continue;
+    }
 
-      mul_v3_v3fl(proxy[vd.i], grab_delta, fade);
+    const float dist = ss->cache->geodesic_dists[vd.index];
+    float fade = bstrength * SCULPT_brush_strength_factor(ss,
+                                                          brush,
+                                                          orig_data.co,
+                                                          dist,
+                                                          orig_data.no,
+                                                          NULL,
+                                                          vd.mask ? *vd.mask : 0.0f,
+                                                          vd.index,
+                                                          thread_id);
 
-      if (vd.mvert) {
-        vd.mvert->flag |= ME_VERT_PBVH_UPDATE;
+    if (grab_silhouette) {
+      float silhouette_test_dir[3];
+      normalize_v3_v3(silhouette_test_dir, grab_delta);
+      if (dot_v3v3(ss->cache->initial_normal, ss->cache->grab_delta_symmetry) < 0.0f) {
+        mul_v3_fl(silhouette_test_dir, -1.0f);
       }
+      float vno[3];
+      normal_short_to_float_v3(vno, orig_data.no);
+      fade *= max_ff(dot_v3v3(vno, silhouette_test_dir), 0.0f);
+    }
+
+    mul_v3_v3fl(proxy[vd.i], grab_delta, fade);
+
+    if (vd.mvert) {
+      vd.mvert->flag |= ME_VERT_PBVH_UPDATE;
     }
   }
   BKE_pbvh_vertex_iter_end;
@@ -4005,6 +4009,11 @@ static void do_grab_brush(Sculpt *sd, Object *ob, PBVHNode **nodes, int totnode)
     sculpt_project_v3_normal_align(ss, ss->cache->normal_weight, grab_delta);
   }
 
+  if (SCULPT_stroke_is_first_brush_step(ss->cache)) {
+    ss->cache->geodesic_dists = SCULPT_geodesic_from_vertex_and_symm(
+        sd, ob, SCULPT_active_vertex_get(ss), ss->cache->initial_radius);
+  }
+
   SculptThreadedTaskData data = {
       .sd = sd,
       .ob = ob,
@@ -6644,6 +6653,7 @@ void SCULPT_cache_free(StrokeCache *cache)
   MEM_SAFE_FREE(cache->detail_directions);
   MEM_SAFE_FREE(cache->prev_displacement);
   MEM_SAFE_FREE(cache->limit_surface_co);
+  MEM_SAFE_FREE(cache->geodesic_dists);
 
   if (cache->pose_ik_chain) {
     SCULPT_pose_ik_chain_free(cache->pose_ik_chain);
@@ -9659,6 +9669,234 @@ static void SCULPT_OT_dyntopo_detail_size_edit(wmOperatorType *ot)
   ot->flag = OPTYPE_REGISTER | OPTYPE_UNDO;
 }
 
+#define SCULPT_GEODESIC_VERTEX_NONE -1
+
+static bool sculpt_geodesic_mesh_test_dist_add(
+    MVert *mvert, const int v0, const int v1, const int v2, float *dists, GSet *initial_vertices)
+{
+  if (BLI_gset_haskey(initial_vertices, POINTER_FROM_INT(v0))) {
+    return false;
+  }
+
+  BLI_assert(dists[v1] != FLT_MAX);
+  if (dists[v0] <= dists[v1]) {
+    return false;
+  }
+
+  float dist0;
+  if (v2 != SCULPT_GEODESIC_VERTEX_NONE) {
+    BLI_assert(dists[v2] != FLT_MAX);
+    if (dists[v0] <= dists[v2]) {
+      return false;
+    }
+    dist0 = geodesic_distance_propagate_across_triangle(
+        mvert[v0].co, mvert[v1].co, mvert[v2].co, dists[v1], dists[v2]);
+  }
+  else {
+    float vec[3];
+    sub_v3_v3v3(vec, mvert[v1].co, mvert[v0].co);
+    dist0 = dists[v1] + len_v3(vec);
+  }
+
+  if (dist0 < dists[v0]) {
+    dists[v0] = dist0;
+    return true;
+  }
+
+  return false;
+}
+
+static float *SCULPT_geodesic_mesh_create(Object *ob,
+                                          GSet *initial_vertices,
+                                          const float limit_radius)
+{
+  SculptSession *ss = ob->sculpt;
+  Mesh *mesh = ob->data;
+
+  const int totvert = SCULPT_vertex_count_get(ss);
+  const int totedge = mesh->totedge;
+
+  const float limit_radius_sq = limit_radius * limit_radius;
+
+  MEdge *edges = mesh->medge;
+  MVert *verts = SCULPT_mesh_deformed_mverts_get(ss);
+
+  float *dists = MEM_malloc_arrayN(totvert, sizeof(float), "distances");
+  BLI_bitmap *edge_tag = BLI_BITMAP_NEW(totedge, "edge tag");
+
+  if (!ss->epmap) {
+    BKE_mesh_edge_poly_map_create(&ss->epmap,
+                                  &ss->epmap_mem,
+                                  mesh->medge,
+                                  mesh->totedge,
+                                  mesh->mpoly,
+                                  mesh->totpoly,
+                                  mesh->mloop,
+                                  mesh->totloop);
+  }
+  if (!ss->vemap) {
+    BKE_mesh_vert_edge_map_create(
+        &ss->vemap, &ss->vemap_mem, mesh->medge, mesh->totvert, mesh->totedge);
+  }
+
+  BLI_LINKSTACK_DECLARE(queue, int);
+  BLI_LINKSTACK_DECLARE(queue_next, int);
+
+  BLI_LINKSTACK_INIT(queue);
+  BLI_LINKSTACK_INIT(queue_next);
+
+  for (int i = 0; i < totvert; i++) {
+    if (BLI_gset_haskey(initial_vertices, POINTER_FROM_INT(i))) {
+      dists[i] = 0.0f;
+    }
+    else {
+      dists[i] = FLT_MAX;
+    }
+  }
+
+  BLI_bitmap *affected_vertex = BLI_BITMAP_NEW(totvert, "affected vertex");
+  GSetIterator gs_iter;
+  GSET_ITER (gs_iter, initial_vertices) {
+    const int v = POINTER_AS_INT(BLI_gsetIterator_getKey(&gs_iter));
+    float *v_co = verts[v].co;
+    for (int i = 0; i < totvert; i++) {
+      if (len_squared_v3v3(v_co, verts[i].co) <= limit_radius_sq) {
+        BLI_BITMAP_ENABLE(affected_vertex, i);
+      }
+    }
+  }
+
+  for (int i = 0; i < totedge; i++) {
+    const int v1 = edges[i].v1;
+    const int v2 = edges[i].v2;
+    if (!affected_vertex[v1] && !affected_vertex[v2]) {
+      continue;
+    }
+    if (dists[v1] != FLT_MAX || dists[v2] != FLT_MAX) {
+      BLI_LINKSTACK_PUSH(queue, i);
+    }
+  }
+
+  do {
+    int e;
+    while (e = BLI_LINKSTACK_POP(queue)) {
+      int v1 = edges[e].v1;
+      int v2 = edges[e].v2;
+
+      if (dists[v1] == FLT_MAX || dists[v2] == FLT_MAX) {
+        if (dists[v1] > dists[v2]) {
+          SWAP(int, v1, v2);
+        }
+        sculpt_geodesic_mesh_test_dist_add(
+                verts, v2, v1, SCULPT_GEODESIC_VERTEX_NONE, dists, initial_vertices);
+      }
+
+      if (ss->epmap[e].count != 0) {
+        for (int pi = 0; pi < ss->epmap[e].count; pi++) {
+          const int poly = ss->epmap[e].indices[pi];
+          if (ss->face_sets[poly] <= 0) {
+            continue;
+          }
+          const MPoly *mpoly = &mesh->mpoly[poly];
+
+          for (int li = 0; li < mpoly->totloop; li++) {
+            const MLoop *mloop = &mesh->mloop[li + mpoly->loopstart];
+            const int v_other = mloop->v;
+            if (ELEM(v_other, v1, v2)) {
+              continue;
+            }
+            if (sculpt_geodesic_mesh_test_dist_add(
+                    verts, v_other, v1, v2, dists, initial_vertices)) {
+              for (int ei = 0; ei < ss->vemap[v_other].count; ei++) {
+                const int e_other = ss->vemap[v_other].indices[ei];
+                int ev_other;
+                if (edges[e_other].v1 == v_other) {
+                  ev_other = edges[e_other].v2;
+                }
+                else {
+                  ev_other = edges[e_other].v1;
+                }
+
+                if (e_other != e && !BLI_BITMAP_TEST(edge_tag, e_other) &&
+                    (ss->epmap[e_other].count == 0 || dists[ev_other] != FLT_MAX)) {
+                  if (BLI_BITMAP_TEST(affected_vertex, v_other) ||
+                      BLI_BITMAP_TEST(affected_vertex, ev_other)) {
+                    BLI_BITMAP_ENABLE(edge_tag, e_other);
+                    BLI_LINKSTACK_PUSH(queue_next, e_other);
+                  }
+                }
+              }
+            }
+          }
+        }
+      }
+    }
+
+    for (LinkNode *lnk = queue_next; lnk; lnk = lnk->next) {
+      const int e = POINTER_AS_INT(lnk->link);
+      BLI_BITMAP_DISABLE(edge_tag, e);
+    }
+
+    BLI_LINKSTACK_SWAP(queue, queue_next);
+
+  } while (BLI_LINKSTACK_SIZE(queue));
+
+  BLI_LINKSTACK_FREE(queue);
+  BLI_LINKSTACK_FREE(queue_next);
+  MEM_SAFE_FREE(edge_tag);
+  MEM_SAFE_FREE(affected_vertex);
+
+  return dists;
+}
+
+float *SCULPT_geodesic_distances_create(Object *ob,
+                                        GSet *initial_vertices,
+                                        const float limit_radius)
+{
+  SculptSession *ss = ob->sculpt;
+  const int totvert = SCULPT_vertex_count_get(ss);
+  switch (BKE_pbvh_type(ss->pbvh)) {
+    case PBVH_FACES:
+      return SCULPT_geodesic_mesh_create(ob, initial_vertices, limit_radius);
+    case PBVH_BMESH:
+    case PBVH_GRIDS:
+      return MEM_calloc_arrayN(totvert, sizeof(float), "empty geodesic distances");
+  }
+  BLI_assert(false);
+  return NULL;
+}
+
+float *SCULPT_geodesic_from_vertex_and_symm(Sculpt *sd,
+                                            Object *ob,
+                                            const int vertex,
+                                            const float limit_radius)
+{
+  SculptSession *ss = ob->sculpt;
+  GSet *initial_vertices = BLI_gset_int_new("initial_vertices");
+
+  const char symm = SCULPT_mesh_symmetry_xyz_get(ob);
+  for (char i = 0; i <= symm; ++i) {
+    if (SCULPT_is_symmetry_iteration_valid(i, symm)) {
+      int v = -1;
+      if (i == 0) {
+        v = vertex;
+      }
+      else {
+        float location[3];
+        flip_v3_v3(location, SCULPT_vertex_co_get(ss, vertex), i);
+        v = SCULPT_nearest_vertex_get(sd, ob, location, FLT_MAX, false);
+      }
+      if (v != -1) {
+        BLI_gset_add(initial_vertices, POINTER_FROM_INT(v));
+      }
+    }
+  }
+
+  float *dists = SCULPT_geodesic_distances_create(ob, initial_vertices, limit_radius);
+  BLI_gset_free(initial_vertices, NULL);
+  return dists;
+}
+
 void ED_operatortypes_sculpt(void)
 {
   WM_operatortype_append(SCULPT_OT_brush_stroke);
diff --git a/source/blender/editors/sculpt_paint/sculpt_intern.h b/source/blender/editors/sculpt_paint/sculpt_intern.h
index d1e17c7e59b..9832525ac73 100644
--- a/source/blender/editors/sculpt_paint/sculpt_intern.h
+++ b/source/blender/editors/sculpt_paint/sculpt_intern.h
@@ -361,6 +361,15 @@ float *SCULPT_boundary_automasking_init(Object *ob,
                                         int propagation_steps,
                                         float *automask_factor);
 
+/* Geodesic distances. */
+float *SCULPT_geodesic_distances_create(struct Object *ob,
+                                        struct GSet *initial_vertices,
+                                        const float limit_radius);
+float *SCULPT_geodesic_from_vertex_and_symm(struct Sculpt *sd,
+                                            struct Object *ob,
+                                            const int vertex,
+                                            const float limit_radius);
+
 /* Filters. */
 void SCULPT_filter_cache_init(struct bContext *C, Object *ob, Sculpt *sd, const int undo_type);
 void SCULPT_filter_cache_free(SculptSession *ss);
@@ -950,6 +959,9 @@ typedef struct StrokeCache {
   bool is_rake_rotation_valid;
   struct SculptRakeData rake_data;
 
+  /* Geodesic distances. */
+  float *geodesic_dists;
+
   /* Face Sets */
   int paint_face_set;