分形抖动Voronoi分区

# 分形抖动Voronoi划分 来源：https://www.boristhebrave.com/2026/08/29/fractal-jittered-voronoi-partitions/ 我发现了一种很酷的Voronoi变体，它能为每个Voronoi单元格生成优美的分形边界，且计算效率很高。这似乎特别适合模拟海岸线——海岸线通常也具有类似的分形特性。 https://www.shadertoy.com/view/sfKSDw 让我详细解释一下。“抖动”Voronoi（https://www.boristhebrave.com/docs/sylves/1/articles/grids/jitteredsquaregrid.html）算法首先创建一个无限网格，然后在每个网格方格内随机选取一个点，称为**生成点**。接着基于这些生成点构建Voronoi图。这与Worley噪声（https://en.wikipedia.org/wiki/Worley_noise）原理相同，区别在于算法输出的不是距离值，而是最近的生成点选择，从而对平面进行划分。它有一个显著优势：易于按点计算。对于给定点*p*，只需查找最近的25个网格单元格，通过哈希伪随机数生成器获取各单元格的生成点，然后找出距离最近的那个即可——无需实际构建Voronoi图，也不必处理平面的无限性问题。 以下是我们的新流程，数学描述如下： 1) 同样从无限方格网格开始，通过PRNG为每个方格选取一个生成点。这些是第0层的根生成点。 2) 创建新网格，方格尺寸减半，为它们选取第1层生成点。对每个第1层生成点，找到其在第0层的父生成点（定义为最近的第0层生成点）。 3) 无限重复步骤2。每次将网格尺寸减半，并在上层中寻找父生成点。 每个生成点所属的划分区域，即为沿父链向上追溯得到的根生成点。随着递归进行，平面会被生成点逐步填充，最终可扩展为平面的几乎处处划分（https://en.wikipedia.org/wiki/Almost_everywhere）。但在实践中如何计算呢？我们既无法处理无限网格，也不能进行无限递归。对于无限网格，可以采用常规抖动Voronoi的处理方式——对于给定点*p*，每个层级中实际相关的单元格数量都是有限的。至于无限递归，显然的解决方案是执行有限次迭代（取决于缩放级别），并构建最终层的Voronoi图，这能很好地近似分形形态。 --- 下面用代码实现。我稍微简化了循环终止条件的写法。 ``` from functools import cache from math import sqrt # 各层级的抖动Voronoi def cell_size(layer): return 2 ** (-layer) @cache def site(layer, cell): """ 返回属于整数网格单元(x, y)的生成点。 hash2必须根据(layer, cell)确定性地返回两个[0, 1)范围内的伪随机数。 """ s = cell_size(layer) offset = hash2(layer, cell) return s * (cell + offset) def nearest_site(layer, p): """查找层级layer中距离点p最近的生成点。""" s = cell_size(layer) centre_cell = floor(p / s) best_cell = None best_site = None best_distance = infinity # 搜索最近的5x5单元格就足够了： # 更远单元格不可能包含最近生成点。 for dx in range(-2, 3): for dy in range(-2, 3): cell = centre_cell + (dx, dy) q = site(layer, cell) d = squared_distance(p, q) if d < best_distance: best_distance = d best_cell = cell best_site = q return best_cell, best_site # 父级/根级关系 @cache def parent(layer, cell): """返回生成点在layer层级的父生成点。""" return nearest_site(layer - 1, site(layer, cell))[0] @cache def root(layer, cell): """返回某层单元格对应的第0层单元格。""" while layer > 0: cell = parent(layer, cell) layer -= 1 return cell def partition(p, depth): """分形划分的有限深度近似计算。""" cell, _ = nearest_site(depth, p) return root(depth, cell) ``` 另一种变体是设置提前终止条件，而非固定深度。如果点*p*恰好非常接近第*i*层的某个生成点，我们可以停止递归并从该处向上追溯父链。这已足够接近生成点，后续所有层级都会认同该结果。 ``` def nearest_site_safe(layer, p, check_safe=False): """ 查找层级layer中距离点p最近的生成点。 同时判断所有可能的更深层生成点是否具有相同的根生成点。 如果是，则进一步递归不会改变结果。 """ s = cell_size(layer) centre_cell = floor(p / s) best_cell = None best_site = None best_distance = infinity # 在任意深层选取的生成点，其第layer层 # 祖先距离p最多为该值。 influence_radius = 2 * sqrt(2) * s possible_roots = set() # ±2范围足以找到最近生成点。 # 检查是否安全停止时使用±3，包含所有可能 # 其生成点位于影响半径内的单元格。 radius = 3 if check_safe else 2 for dx in range(-radius, radius + 1): for dy in range(-radius, radius + 1): cell = centre_cell + (dx, dy) q = site(layer, cell) d = squared_distance(p, q) if d < best_distance: best_distance = d best_cell = cell best_site = q if check_safe and d <= influence_radius ** 2: possible_roots.add(root(layer, cell)) safe = len(possible_roots) == 1 return best_cell, best_site, safe def partition_adaptive(p, max_depth): """ 计算包含点p的划分区域，当保证进一步递归不会改变结果时提前终止。 若无法达成此保证，则使用max_depth作为有限近似。 """ for layer in range(max_depth + 1): cell, _, safe = nearest_site(layer, p, check_safe=True) if safe: break return root(layer, cell) ```