如何在不丢弃的情况下在八面体中生成随机点？答案

【问题标题】：How to generate random points within an octahedron without discarding?如何在不丢弃的情况下在八面体中生成随机点？
【发布时间】：2020-03-03 11:31:28
【问题描述】：

我需要一个八面体中的随机点，均匀分布。我将八面体定义为所有点都满足abs(x) + abs(y) + abs(z) <= 1 的体积，其中 abs 给出绝对值。 IE：六个顶点中的每一个都在一个轴上，距离 0,0,0 1。也许你可以称它为单位八面体。

考虑到定义，我可以天真地生成一个点，如下所示：

val x: Double = nextDouble() // 0-1 range
val y = nextDouble(1.0 -x) // 1-x is upper bound, probably <1
val z = nextDouble(1.0 -(x+y))

问题在于这倾向于较小的 y 值和较小的 z 值。显然不是均匀分布。同样清楚的是，所有这些点都在八个象限之一中。

我避免使用 discard 方法，因为这个函数会被大量调用，而且看起来我应该能够做得比丢掉大部分分数更好。

请注意，dual of the octahedron 是立方体。正因为如此，我有一个暗示，可能存在一个简单的函数来将立方体内的任何点转换为八面体内，但这只是我仍在探索的直觉。

【问题讨论】：

标签： math random 3d polyhedra

【解决方案1】：

这里的想法是——在 D+1 中从Dirichlet distribution 中采样点，选择 D 点，这样你在单纯形中是一致的

x₀+x₁+x₂i >= 0

要制作八面体，请随机选择八分圆来表达您的观点。

Python 代码

import math
import random

def Dirichlet():
    """sample 4d Dirichlet"""
    x0 = -math.log(1.0 - random.random()) # exponential
    x1 = -math.log(1.0 - random.random()) # exponential
    x2 = -math.log(1.0 - random.random()) # exponential
    x3 = -math.log(1.0 - random.random()) # exponential
    s = 1.0/(x0+x1+x2+x3) # scaling

    return (x0*s, x1*s, x2*s, x3*s)

def OctahedronSampling():

    x0, x1, x2, _ = Dirichlet()

    octant = random.randint(0, 7)

    if octant == 0:
        return (x0, x1, x2)
    elif octant == 1:
        return (x0, -x1, x2)
    elif octant == 2:
        return (x0, x1, -x2)
    elif octant == 3:
        return (x0, -x1, -x2)
    elif octant == 4:
        return (-x0, x1, x2)
    elif octant == 5:
        return (-x0, -x1, x2)
    elif octant == 6:
        return (-x0, x1, -x2)
    elif octant == 7:
        return (-x0, -x1, -x2)

    return None

for _ in range(0, 2000):
    x0, x1, x2 = OctahedronSampling()

    print(f"{x0}   {x1}   {x2}")

这是 2K 点的快速图表

【讨论】：

【解决方案2】：

你知道如何在一个均匀分布的立方体中选择点，一个立方体可以被分解成八个方形金字塔。（对不起，我无法提供图形。）

我会从一个立方体开始：abs(x) <= 1; abs(y) <= 1; abs(z) <= 1

在其中选择一个点（列向量，(x, y, z)），然后反射将其带入“顶部和底部”金字塔：

if abs(x) > abs(z), swap x and z. Equivalently, multiply by

0 0 1
0 1 0
1 0 0

if abs(y) > abs(z), swap y and z. Equivalently, multiply by 

1 0 0
0 0 1
0 1 0

然后将两个金字塔倒转成一个八面体：

if z>0
 z = 1-z

if z<0
 z = -1-z

然后旋转和缩放：

multiply by

1/2 -1/2  0
1/2  1/2  0
  0    0  1

【讨论】：