生成随机数a,b,c,使得a ^ 2 + b ^ 2 + c ^ 2 = 1

fou*_*nes 6 python random

要在Python中进行一些模拟,我正在尝试生成数字a,b,c a^2 + b^2 + c^2 = 1.我认为a在0和1之间生成一些,然后b在0和之间生成一些sqrt(1 - a^2),然后c= sqrt(1 - a^2 - b^2)将起作用.

浮点值很好,平方和应该接近 1.我想在一些迭代中继续生成它们.

作为Python的新手,我不确定如何做到这一点.允许否定.

编辑:非常感谢您的答案!

Joh*_*anL 8

根据stats.stackexchange.com 上的这个答案,您应该使用正态分布值来在球体上获得均匀分布的值.这意味着,你可以这样做:

import numpy as np
abc = np.random.normal(size=3)
a,b,c = abc/np.sqrt(sum(abc**2))
Run Code Online (Sandbox Code Playgroud)


MSe*_*ert 8

如果您对概率密度感兴趣,我决定对不同方法进行比较:

import numpy as np
import random
import math

def MSeifert():
    a = 1
    b = 1
    while a**2 + b**2 > 1:  # discard any a and b whose sum of squares already exceeds 1
        a = random.random()
        b = random.random()
    c = math.sqrt(1 - a**2 - b**2)  # fixed c
    return a, b, c

def VBB():
    x = np.random.uniform(0,1,3) # random numbers in [0, 1)
    x /= np.sqrt(x[0] ** 2 + x[1] ** 2 + x[2] ** 2)
    return x[0], x[1], x[2]

def user3684792():
    theta = random.uniform(0, 0.5*np.pi)
    phi = random.uniform(0, 0.5*np.pi)
    return np.sin(theta)* np.cos(phi), np.sin(theta)*np.sin(phi), np.cos(theta)

def JohanL():
    abc = np.random.normal(size=3)
    a,b,c = abc/np.sqrt(sum(abc**2))
    return a, b, c

def SeverinPappadeux():
    cos_th = 2.0*random.uniform(0, 1.0) - 1.0
    sin_th = math.sqrt(1.0 - cos_th*cos_th)
    phi = random.uniform(0, 2.0*math.pi)
    return sin_th * math.cos(phi), sin_th * math.sin(phi), cos_th
Run Code Online (Sandbox Code Playgroud)

并绘制分布:

%matplotlib notebook

import matplotlib.pyplot as plt

f, axes = plt.subplots(3, 4)

for func_idx, func in enumerate([MSeifert, JohanL, user3684792, VBB]):
    axes[0, func_idx].set_title(str(func.__name__))
    res = [func() for _ in range(50000)]
    for idx in range(3):
        axes[idx, func_idx].hist([i[idx] for i in res], bins='auto')

axes[0, 0].set_ylabel('a')
axes[1, 0].set_ylabel('b')
axes[2, 0].set_ylabel('c')

plt.tight_layout()
Run Code Online (Sandbox Code Playgroud)

结果如下:

在此输入图像描述

说明:行显示分布a,b并且c列分别显示不同方法的直方图(分布).

在范围内给出均匀随机分布的唯一方法(-1, 1)是JohanLs和Severin Pappadeux的方法.所有其他方法都有一些功能,如尖峰或范围内的功能行为[0, 1).请注意,这两个解决方案目前的值介于-1和1之间,而所有其他方法的值介于0和1之间.

  • 好地块.值得指出的是,球体上的均匀分布给出均匀的边缘分布这一事实相当微妙,并且对于2球体是特殊的:相应的结果不适用于圆形或超球形. (2认同)

use*_*792 6

我认为它实际上是一个很酷的问题,一个很好的方法就是使用球形极坐标并随机生成角度.

import random
import numpy as np
def random_pt():
    theta = random.uniform(0, 0.5*np.pi)
    phi = random.uniform(0, 0.5*np.pi)
    return np.sin(theta)* np.cos(phi), np.sin(theta)*np.sin(phi), np.cos(theta)
Run Code Online (Sandbox Code Playgroud)

  • 这里有一个问题:点将在极点附近聚集. (2认同)

VBB*_*VBB 3

“正确”的答案取决于您是否正在寻找空间中、球体表面或其他物体上的均匀随机分布。如果您正在寻找球体表面上的点,您仍然需要担心cos(theta)会导致点在球体两极附近出现“聚集”的因素。由于您的问题尚不清楚确切的性质,因此这里有一个应该有效的“完全随机”分布:

x = np.random.uniform(0,1,3) # random numbers in [0, 1)
x /= np.sqrt(x[0] ** 2 + x[1] ** 2 + x[2] ** 2)
Run Code Online (Sandbox Code Playgroud)

x = np.random.uniform(0, 1, (3, n))这里的另一个优点是,由于我们使用 numpy 数组,因此您也可以通过使用for any来快速扩展到大型点集n