【2024认证杯小美赛B题】完整解析与代码分享（全网最全）

（多重覆盖）：一个大的有限赫希数

背景

被称为平面镶嵌的问题涉及使用单个多边形,该多边形经过旋转、平移和镜像处理,以在平面上实现无间隙和无重叠的镶嵌。有些多边形能够实现镶嵌,例如正方形,而其他多边形则不能,例如五边形。我们现在将重点关注那些无法实现镶嵌的多边形。
举例来说,当一个五边形位于一个平面内时,使用同一个五边形来覆盖其整个周长而不留下间隙是不可能的。然而,有些多边形可以通过某种方式构建,使其覆盖周长(日冤)而不留下间隙,但不能继续。海因里希·赫希在1968年构建的图形就是一个例子。

此外,存在一些多边形,它们可以扩展以涵盖额外一周的边缘,但无法进一步扩展。这可以通过安妮·方丹在1991年构建的图形来说明[1]。

在这种情况下,能够填满的最大圆数被称为赫希数。
正式的数学定义如下[2]:平面的镶嵌被定义为将平面划分为较小的区域或瓷砖。瓷砖的第零个日冤被定义为瓷砖本身。对于k0,第k个日冤被定义为与第(k-1)个日冤共享边界点的瓷砖集合。给定图形S的海斯数被定义为存在平面到较小区域(称为瓷砖)的镶嵌的最大k值。图形S的海斯数被定义为存在平面的镶嵌以及该镶嵌中的瓷砖t,使得t的第零到第k个日冤中的所有瓷砖都是简单连通区域的最大k值。在某些关于这个问题的研究中,这个定义被修改,要求瓷砖的第零到第k个日冤的并集是一个简单连通区域。

赫希问题的目标是确定理论上的最大赫希数。已知的最大学赫希数为6,
如下图所示,这是博扬·巴希奇在2020年构建的[3]。

任务

1. 目标是构建一个高效的数学模型和算法,能够生成具有最大可能赫希
   数的多边形。本文应详细描述算法的原理和实现过程,以及最终结果。
   不必达到或超过已知的最优结果;然而,如果算法能自动生成最大可
   能的结果,则令人满意。
2. 请对该算法的复杂度上限进行估算。

平面镶嵌理论中的 Heesch 数

原理

平面镶嵌指的是用一种形状（比如正方形、三角形、五边形）去填满整个平面，就像拼图一样，不留任何缝隙，也不重叠。这在日常生活中很常见，比如瓷砖铺地、蜂巢结构。

但有些形状并不能填满整个平面，比如一个普通的五边形。虽然这些形状不能无限延伸镶嵌，但有些是可以围绕自己“镶嵌几圈”的。


具体来说：

1. 第一圈：用多个一样的形状围绕一个形状，像把它包裹起来。
2. 第二圈：在第一圈的外侧继续包裹一圈。
3. 第三圈：继续在外面再镶嵌一圈……
   这样一圈一圈地包裹，直到某一圈没办法继续镶嵌下去为止。这个“最大能包裹的圈数”就被称为 Heesch 数。

举例：

• 正方形：正方形可以一直无缝镶嵌下去，它的 Heesch 数是“无限”。
• 普通五边形：某些五边形可能只能围绕自己包裹 2 圈，第三圈开始就会出现缝隙或重叠，那么它的 Heesch 数是 2。
• 复杂的多边形：有些形状可能能包裹 5 或 6 圈，超过了以后就无法继续了。

算法步骤

Step 1: 初始多边形生成

• 生成规则多边形（如六边形）作为起始形状。
• 引入形状变形（通过顶点位置的扰动）生成候选多边形。

Step 2: 镶嵌模拟

• 镶嵌生成：
  • 模拟第 0 到第 k 环镶嵌，通过平移、旋转、镜像生成下一环。
• 连通性检测：
  • 使用几何库（如 Shapely）计算镶嵌区域是否简单连通。
• 冲突检测：
  • 检查多边形之间是否存在重叠或缝隙。

Step 3: 形状优化

• 优化目标：最大化 Heesch 数。
• 优化方法：
  • 使用粒子群算法优化顶点坐标 V。
  • 目标函数：f(V)=k，其中 k 是当前多边形的 Heesch 数。
• 搜索策略：
  • 形状调整后，重新评估 Heesch 数，保留改进结果。

Step 4: 结果输出

• 输出最终的多边形形状和对应的最大 Heesch 数。

粒子群算法

代码实现分析

import matplotlib.pyplot as plt
from shapely.geometry import Polygon, MultiPolygon
from shapely.affinity import rotate, translate
import numpy as np
import random

def create_polygon(vertices):
    """
    根据顶点创建多边形。
    """
    return Polygon(vertices)

def generate_tiling(polygon, layers):
    """
    模拟多边形镶嵌。
    """
    current_layer = [polygon]
    all_layers = [polygon]
    
    for _ in range(layers):
        next_layer = []
        for poly in current_layer:
            for angle in range(0, 360, 360 // len(polygon.exterior.coords) - 1):
                rotated = rotate(polygon, angle, origin=(0, 0))
                for x, y in poly.exterior.coords:
                    translated = translate(rotated, xoff=x, yoff=y)
                    if not any(p.intersects(translated) for p in all_layers):
                        next_layer.append(translated)
        if not next_layer:
            break
        current_layer = next_layer
        all_layers += next_layer
    return MultiPolygon(all_layers)

def is_connected(multipolygon):
    """
    检查镶嵌区域是否简单连通。
    """
    union = multipolygon.buffer(0.01)
    return union.is_valid and union.is_simple

def compute_heesch_number(vertices, max_layers=10):
    """
    计算给定多边形顶点的 Heesch 数。
    """
    polygon = create_polygon(vertices)
    tiling = generate_tiling(polygon, max_layers)
    if is_connected(tiling):
        return len(tiling.geoms) // len(polygon.exterior.coords)
    return 0

def particle_swarm_optimization(num_particles, num_iterations, vertices_dim, max_layers):
    """
    基于粒子群优化的多边形优化算法。
    """
    # 粒子初始化
    particles = [np.random.uniform(-1, 1, (vertices_dim, 2)) for _ in range(num_particles)]
    velocities = [np.zeros((vertices_dim, 2)) for _ in range(num_particles)]
    personal_best = particles[:]
    personal_best_scores = [compute_heesch_number(p, max_layers) for p in particles]
    global_best = personal_best[np.argmax(personal_best_scores)]
    global_best_score = max(personal_best_scores)

    # 超参数
    inertia = 0.5  # 惯性权重
    cognitive = 1.5  # 自我认知权重
    social = 1.5  # 社会认知权重

    for iteration in range(num_iterations):
        for i in range(num_particles):
            # 更新粒子速度和位置
            r1, r2 = np.random.random(), np.random.random()
            velocities[i] = (
                inertia * velocities[i]
                + cognitive * r1 * (personal_best[i] - particles[i])
                + social * r2 * (global_best - particles[i])
            )
            particles[i] += velocities[i]

            # 计算适应度
            score = compute_heesch_number(particles[i], max_layers)
            if score > personal_best_scores[i]:
                personal_best[i] = particles[i]
                personal_best_scores[i] = score

            if score > global_best_score:
                global_best = particles[i]
                global_best_score = score

        print(f"Iteration {iteration + 1}/{num_iterations}, Best Heesch Number: {global_best_score}")

    return global_best, global_best_score

def visualize_polygon(vertices, heesch_number):
    """
    可视化最终多边形及其镶嵌。
    """
    polygon = create_polygon(vertices)
    tiling = generate_tiling(polygon, heesch_number)
    
    plt.figure(figsize=(8, 8))
    plt.gca().set_aspect('equal', adjustable='box')
    
    for poly in tiling.geoms:
        x, y = poly.exterior.xy
        plt.fill(x, y, alpha=0.6)
    plt.show()

# 初始多边形顶点维度
vertices_dim = 6  # 六边形顶点
num_particles = 30  # 粒子数量
num_iterations = 100  # 迭代次数
max_layers = 10  # 最大镶嵌层数

# 粒子群优化
best_vertices, best_heesch_number = particle_swarm_optimization(num_particles, num_iterations, vertices_dim, max_layers)
print(f"最大 Heesch 数: {best_heesch_number}")

# 可视化结果
visualize_polygon(best_vertices, best_heesch_number)

完整代码

请看下方~