如何生成一个随机的凸多边形?

17

我正在尝试设计一种生成随机二维凸多边形的方法。它必须具备以下特点:

  • 坐标应为整数;
  • 多边形应位于以(0,0)和(C,C)为角的正方形内,其中C已知;
  • 多边形的顶点数应接近给定的数字N。

例如,生成具有10个顶点且位于区间[0..100]x[0..100]的随机多边形。

这项任务的难点在于坐标必须为整数。

我的方法是在给定的正方形中生成一组随机点,并计算这些点的凸包。但是得到的凸包顶点很少,远少于N。

有什么想法吗?

9个回答

8
这是我所知道的最快算法,可以以相等的概率生成每个凸多边形。输出正好有N个顶点,运行时间为O(N log N),因此即使是大型多边形也可以快速生成。
生成两个列表 X 和 Y,包含 N 个介于 0 和 C 之间的随机整数,确保没有重复项。 对 X 和 Y 进行排序并存储它们的最大值和最小值。 将其他元素(不是最大值或最小值)随机分成两组:X1 和 X2,以及 Y1 和 Y2。 在这些列表的开头和结尾重新插入最小值和最大值元素(minX 在 X1 和 X2 的开头,maxX 在末尾等)。 找到连续的差异(X1[i + 1] - X1[i]),对于第二组反转顺序(X2[i] - X2[i + 1])。将它们存储在列表 XVec 和 YVec 中。 随机打乱 YVec,并将每对 XVec[i] 和 YVec[i] 视为 2D 向量。 按角度对这些向量进行排序,然后将它们端对端排列形成一个多边形。 将多边形移动到原始的最小和最大坐标处。

这里提供了一个动画和Java实现: 生成随机凸多边形

该算法基于Pavel Valtr的一篇论文:“n个随机点在凸位置的概率。”《离散与计算几何》13.1 (1995): 637-643。


如果N == 3,那么x1或x2将必须为空,你应该怎么做? - Andrew
1
X1和X2在下一步后都包括最小值和最大值元素。中间元素随机添加到其中一个。您基本上正在决定中点是三角形顶部还是底部。当然,对于N=3来说这些都是多余的,因为任何三个点都形成唯一的三角形,且所有三角形都是凸的。 - Mangara

4

根据@Mangara的答案,这里有一个JAVA实现,如果有人对Python移植感兴趣。

import random
from math import atan2


def to_convex_contour(vertices_count,
                      x_generator=random.random,
                      y_generator=random.random):
    """
    Port of Valtr algorithm by Sander Verdonschot.

    Reference:
        http://cglab.ca/~sander/misc/ConvexGeneration/ValtrAlgorithm.java

    >>> contour = to_convex_contour(20)
    >>> len(contour) == 20
    True
    """
    xs = [x_generator() for _ in range(vertices_count)]
    ys = [y_generator() for _ in range(vertices_count)]
    xs = sorted(xs)
    ys = sorted(ys)
    min_x, *xs, max_x = xs
    min_y, *ys, max_y = ys
    vectors_xs = _to_vectors_coordinates(xs, min_x, max_x)
    vectors_ys = _to_vectors_coordinates(ys, min_y, max_y)
    random.shuffle(vectors_ys)

    def to_vector_angle(vector):
        x, y = vector
        return atan2(y, x)

    vectors = sorted(zip(vectors_xs, vectors_ys),
                     key=to_vector_angle)
    point_x = point_y = 0
    min_polygon_x = min_polygon_y = 0
    points = []
    for vector_x, vector_y in vectors:
        points.append((point_x, point_y))
        point_x += vector_x
        point_y += vector_y
        min_polygon_x = min(min_polygon_x, point_x)
        min_polygon_y = min(min_polygon_y, point_y)
    shift_x, shift_y = min_x - min_polygon_x, min_y - min_polygon_y
    return [(point_x + shift_x, point_y + shift_y)
            for point_x, point_y in points]


def _to_vectors_coordinates(coordinates, min_coordinate, max_coordinate):
    last_min = last_max = min_coordinate
    result = []
    for coordinate in coordinates:
        if _to_random_boolean():
            result.append(coordinate - last_min)
            last_min = coordinate
        else:
            result.append(last_max - coordinate)
            last_max = coordinate
    result.extend((max_coordinate - last_min,
                   last_max - max_coordinate))
    return result


def _to_random_boolean():
    return random.getrandbits(1)

2

以下翻译仅供参考,可能不完全,请您谅解。

如果N<3,则退出。生成一个具有N个顶点的单位圆,并将其随机旋转[0..90]度。

随机从原点向外挤出每个顶点,并使用相邻顶点与原点之间的叉积符号来确定凸性。这是速度和质量之间存在折衷的步骤。

设置顶点后,找到距离原点最远的顶点。将每个顶点除以该大小以归一化多边形,然后按比例缩放(C/2)。平移至(C/2,C/2)并转换为整数。


这里有一个凸性测试的更多信息链接: http://www.gamedev.net/topic/561441-polygon-convexity-test-cant-get-it-right/ - scgrn
这适用于浮点坐标。你如何确保转换回整数时,多边形不会变成凹面的?你可以删除有问题的顶点,但这可能会大幅减少顶点数量。 - Jasiu

1
这是使用numpy的Valtr算法的另一版本。 :)
import numpy as np
import numpy.typing and npt


def generateConvex(n: int) -> npt.NDArray[np.float64]:
    '''
    Generate convex shappes according to Pavel Valtr's 1995 alogrithm. Ported from
    Sander Verdonschot's Java version, found here:
    https://cglab.ca/~sander/misc/ConvexGeneration/ValtrAlgorithm.java
    '''
    # initialise random coordinates
    X_rand, Y_rand = np.sort(np.random.random(n)), np.sort(np.random.random(n))
    X_new, Y_new = np.zeros(n), np.zeros(n)

    # divide the interior points into two chains
    last_true = last_false = 0
    for i in range(1, n):
        if i != n - 1:
            if random.getrandbits(1):
                X_new[i] = X_rand[i] - X_rand[last_true]
                Y_new[i] = Y_rand[i] - Y_rand[last_true]
                last_true = i
            else:
                X_new[i] = X_rand[last_false] - X_rand[i]
                Y_new[i] = Y_rand[last_false] - Y_rand[i]
                last_false = i
        else:
            X_new[0] = X_rand[i] - X_rand[last_true]
            Y_new[0] = Y_rand[i] - Y_rand[last_true]
            X_new[i] = X_rand[last_false] - X_rand[i]
            Y_new[i] = Y_rand[last_false] - Y_rand[i]

    # randomly combine x and y and sort by polar angle
    np.random.shuffle(Y_new)
    vertices = np.stack((X_new, Y_new), axis=-1)
    vertices = vertices[np.argsort(np.arctan2(vertices[:, 1], vertices[:, 0]))]

    # arrange points end to end to form a polygon
    vertices = np.cumsum(vertices, axis=0)

    # center around the origin
    x_max, y_max = np.max(vertices[:, 0]), np.max(vertices[:, 1])
    vertices[:, 0] += ((x_max - np.min(vertices[:, 0])) / 2) - x_max
    vertices[:, 1] += ((y_max - np.min(vertices[:, 1])) / 2) - y_max

    return vertices

1
一个简单的算法如下:
  1. 从随机线开始(一个由两个顶点和两条边组成的多边形)
  2. 选择多边形的一条随机边E
  3. 在这条边上生成一个新的随机点P
  4. 取一条垂直于E且经过点P的线L。通过计算线T与相邻于E的两条边所定义的线之间的交点,计算当凸度不被破坏时P的最大偏移量。
  5. 在该范围内随机偏移点P。
  6. 如果点数不足,则从第2步重新开始。

那怎么支持整数坐标? - Jasiu
@Jasiu:我不明白为什么它不支持整数坐标。只需将所有生成的点计算为整数,然后将结果夹紧到所需范围即可。 - Juraj Blaho

1
这里是C++11版本的Pavel Valtr算法实现,该算法源自于Mangara的回答,其中加入了一些类似lewiswolf的回答的技巧,并通过X和Y坐标的分离处理增加了更多的随机性。请注意,保留了HTML标记。
#include <algorithm>
#include <iostream>
#include <random>

struct randPoly {
  int RND_MAX = 655369;
  std::random_device dev;
  std::mt19937 rng;
  std::uniform_int_distribution<std::mt19937::result_type> random_numer;
  std::uniform_int_distribution<std::mt19937::result_type> random_logic;

  randPoly() : rng(dev()), random_numer(0, RND_MAX), random_logic(0, 1) {}
  virtual ~randPoly() {}

  int generate(const int n, const double r0, std::vector<double>& poly_x,
               std::vector<double>& poly_y) {
    auto gen = [&]() { return random_numer(rng); };
    // initialize random samples and sort them
    int m = n / 2;
    std::vector<int> x(n), y(n), vx(n), vy(n), idx(n);
    std::vector<double> a(n);
    std::generate(x.begin(), x.end(), gen);
    std::generate(y.begin(), y.end(), gen);
    std::iota(idx.begin(), idx.end(), 0);
    std::sort(x.begin(), x.end());
    std::sort(y.begin(), y.end());
    // divide samples and get vector component
    int x0 = x[0], x1 = x0;
    for (int k = 1; k < n - 1; ++k) {
      if (random_logic(rng)) {
        vx[k - 1] = x[k] - x0;
        x0 = x[k];
      } else {
        vx[k - 1] = x1 - x[k];
        x1 = x[k];
      }
    }
    vx[n - 2] = x[n - 1] - x0;
    vx[n - 1] = x1 - x[n - 1];
    int y0 = y[0], y1 = y0;
    for (int k = 1; k < n - 1; ++k) {
      if (random_logic(rng)) {
        vy[k - 1] = y[k] - y0;
        y0 = y[k];
      } else {
        vy[k - 1] = y1 - y[k];
        y1 = y[k];
      }
    }
    vy[n - 2] = y[n - 1] - y0;
    vy[n - 1] = y1 - y[n - 1];
    // random pair up vector components and sort by angle
    std::shuffle(vy.begin(), vy.end(), rng);
    for (int k = 0; k < n; ++k) {
      a[k] = std::atan2(vy[k], vx[k]);
    }
    std::sort(idx.begin(), idx.end(),
              [&a](int& lhs, int& rhs) { return a[lhs] < a[rhs]; });
    // form the polygon by connencting vectors
    double x_max = 0, y_max = 0, x_min = 0, y_min = 0;
    x[0] = y[0] = 0;
    for (int k = 1; k < n; ++k) {
      x[k] = x[k - 1] + vx[idx[k - 1]];
      y[k] = y[k - 1] + vy[idx[k - 1]];
      if (x[k] > x_max) {
        x_max = x[k];
      } else if (x[k] < x_min) {
        x_min = x[k];
      }
      if (y[k] > y_max) {
        y_max = y[k];
      } else if (y[k] < y_min) {
        y_min = y[k];
      }
    }
    // center and resize the polygon
    poly_x.resize(n);
    poly_y.resize(n);
    double x_offset = -(x_max + x_min) / 2.0;
    double y_offset = -(y_max + y_min) / 2.0;
    double scale = r0 / std::max(x_max - x_min, y_max - y_min);
    for (int k = 0; k < n; ++k) {
      poly_x[k] = scale * (x[k] + x_offset);
      poly_y[k] = scale * (y[k] + y_offset);
    }
    return 0;
  }
};

int main(int, char**) {
  randPoly rp;
  std::vector<double> poly_x, poly_y;
  rp.generate(8, 2.0, poly_x, poly_y);
  for (int k = 0; k < poly_x.size(); ++k) {
    std::cout << poly_x[k] << "  " << poly_y[k] << std::endl;
  }
}

Example shown in Rviz


1

多亏了@Mangara的回答@Azat的回答,我也制作了Ruby端口:

#!/usr/bin/env ruby
# frozen_string_literal: true

module ValtrAlgorithm
  module_function def random_polygon(length)
    raise ArgumentError, "length should be > 2" unless length > 2

    min_x, *xs, max_x = Array.new(length) { rand }.sort
    min_y, *ys, max_y = Array.new(length) { rand }.sort
    # Divide the interior points into two chains and
    # extract the vector components.
    vec_xs = to_random_vectors(xs, min_x, max_x)
    vec_ys = to_random_vectors(ys, min_y, max_y).
      # Randomly pair up the X- and Y-components
      shuffle
    # Combine the paired up components into vectors
    vecs = vec_xs.zip(vec_ys).
      # Sort the vectors by angle, in a counter clockwise fashion. Remove the
      # `-` to make it clockwise.
      sort_by { |x, y| - Math.atan2(y, x) }

    # Lay them end-to-end
    point_x = point_y = 0
    min_polygon_x = min_polygon_y = 0
    points = []
    vecs.each do |vec_x, vec_y|
      points.append([vec_x, vec_y])
      point_x += vec_x
      point_y += vec_y
      min_polygon_x = [min_polygon_x, point_x].min
      min_polygon_y = [min_polygon_y, point_y].min
    end
    shift_x = min_x - min_polygon_x
    shift_y = min_y - min_polygon_y
    result = points.map { |point_x, point_y| [point_x + shift_x, point_y + shift_y] }
    # Append first point to make it a valid linear ring
    result << result.first
  end

  private def to_random_vectors(coordinates, min, max)
    last_min = last_max = min
    ary = []
    coordinates.each do |coordinate|
      if rand > 0.5
        ary << coordinate - last_min
        last_min = coordinate
      else
        ary << last_max - coordinate
        last_max = coordinate
      end
    end
    ary << max - last_min << last_max - max
  end
end

0

你的初始方法是正确的 - 计算凸包是满足随机性、凸性和整数性的唯一方法。

我能想到的优化算法以获得更多点的唯一方法是将它们围绕一个圆形排列,而不是完全随机。你的点应该更可能靠近正方形的“边缘”而不是中心。在中心,概率应该是 ~0,因为多边形必须是凸的。

一个简单的选择是设置点出现的最小半径 - 可能是C/2或C*0.75。计算C正方形的中心,如果一个点太靠近,就将其移离中心直到达到最小距离。


0
我认为上面使用Numpy的建议非常好,但是可以通过避免使用for循环来改进generateConvex()函数。 另外,我想添加一些Matplotlib的代码来可视化结果。
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.patches import Polygon
rng = np.random.default_rng()

def generateConvex(n):    
    # initialise random coordinates
    XY_rand = np.sort(rng.random((n, 2)), axis=0)
    
    # divide the interior points into two chains
    rand_bool = rng.choice([True, False], n-2)
    pos, neg = XY_rand[1:-1][rand_bool], XY_rand[1:-1][~rand_bool]
    
    pos = np.vstack((XY_rand[0], pos, XY_rand[-1]))
    neg = np.vstack((XY_rand[0], neg, XY_rand[-1]))
    vertices = np.vstack((pos[1:] - pos[:-1], neg[:-1] - neg[1:]))
     
    # randomly combine x and y and sort by polar angle
    rng.shuffle(vertices[:,1])
    vertices = vertices[np.argsort(np.arctan2(vertices[:, 1], vertices[:, 0]))]
    
    # arrange points end to end to form a polygon
    vertices = np.cumsum(vertices, axis=0)
    
    # center around the origin
    x_max, y_max = np.max(vertices[:, 0]), np.max(vertices[:, 1])
    vertices[:, 0] -= (x_max + np.min(vertices[:, 0])) / 2
    vertices[:, 1] -= (y_max + np.min(vertices[:, 1])) / 2
    
    return vertices

if __name__ == '__main__':
    n = 42
    p = Polygon(generateConvex(n))
    fig, ax = plt.subplots()

    ax.add_patch(p)
    plt.title(f'{n}-sided convex polygon')
    plt.axis('equal')
    plt.show()

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接