【问题标题】:How can I set a minimum distance constraint for generating points with numpy.random.rand?如何设置最小距离约束以使用 numpy.random.rand 生成点?
【发布时间】:2015-02-14 10:33:37
【问题描述】:

我正在尝试生成一个有效的代码来生成一些随机位置向量,然后我用它们来计算对相关函数。我想知道是否有直接的方法来限制放置在我的盒子中的任何两点之间允许的最小距离。

我目前的代码如下:

def pointRun(number, dr):
"""
Compute the 3D pair correlation function
for a random distribution of 'number' particles
placed into a 1.0x1.0x1.0 box.
"""
## Create array of distances over which to calculate.   
    r = np.arange(0., 1.0+dr, dr)

## Generate list of arrays to define the positions of all points,
##    and calculate number density.
    a = np.random.rand(number, 3)
    numberDensity = len(a)/1.0**3

## Find reference points within desired region to avoid edge effects. 
    b = [s for s in a if all(s > 0.4) and all(s < 0.6) ]

## Compute pairwise correlation for each reference particle
    dist = scipy.spatial.distance.cdist(a, b, 'euclidean')
    allDists = dist[(dist < np.sqrt(3))]

## Create histogram to generate radial distribution function, (RDF) or R(r)
    Rr, bins = np.histogram(allDists, bins=r, density=False)

## Make empty containers to hold radii and pair density values.
    radii = []
    rhor = []

## Normalize RDF values by distance and shell volume to get pair density.
    for i in range(len(Rr)):
        y = (r[i] + r[i+1])/2.
        radii.append(y)
        x = np.average(Rr[i])/(4./3.*np.pi*(r[i+1]**3 - r[i]**3))
        rhor.append(x)

## Generate normalized pair density function, by total number density
    gr = np.divide(rhor, numberDensity)
    return radii, gr

我之前尝试过使用一个循环来计算每个点的所有距离,然后接受或拒绝。如果我使用很多点,这种方法很慢。

【问题讨论】:

  • 您也许可以使用this algorithm 的变体,但需要注意的是它具有最低分辨率(分辨率越高,成本越高)。

标签: python numpy random distribution correlation


【解决方案1】:

基于@Samir 的回答,为了您的方便,将其设为可调用函数:)

import numpy as np
import matplotlib.pyplot as plt

def generate_points_with_min_distance(n, shape, min_dist):
    # compute grid shape based on number of points
    width_ratio = shape[1] / shape[0]
    num_y = np.int32(np.sqrt(n / width_ratio)) + 1
    num_x = np.int32(n / num_y) + 1

    # create regularly spaced neurons
    x = np.linspace(0., shape[1]-1, num_x, dtype=np.float32)
    y = np.linspace(0., shape[0]-1, num_y, dtype=np.float32)
    coords = np.stack(np.meshgrid(x, y), -1).reshape(-1,2)

    # compute spacing
    init_dist = np.min((x[1]-x[0], y[1]-y[0]))

    # perturb points
    max_movement = (init_dist - min_dist)/2
    noise = np.random.uniform(low=-max_movement,
                                high=max_movement,
                                size=(len(coords), 2))
    coords += noise

    return coords

coords = generate_points_with_min_distance(n=8, shape=(2448,2448), min_dist=256)

# plot
plt.figure(figsize=(10,10))
plt.scatter(coords[:,0], coords[:,1], s=3)
plt.show()

【讨论】:

    【解决方案2】:

    这是一个使用 numpy 的可扩展 O(n) 解决方案。它的工作原理是首先指定一个等距的点网格,然后对这些点进行一定程度的扰动,使点之间的距离最多保持min_dist

    您需要调整点数、盒子形状和扰动灵敏度以获得所需的min_dist

    注意:如果您固定一个框的大小并指定每个点之间的最小距离,那么您可以绘制的满足最小距离的点的数量是有限制的。 .

    import numpy as np
    import matplotlib.pyplot as plt
    
    # specify params
    n = 500
    shape = np.array([64, 64])
    sensitivity = 0.8 # 0 means no movement, 1 means max distance is init_dist
    
    # compute grid shape based on number of points
    width_ratio = shape[1] / shape[0]
    num_y = np.int32(np.sqrt(n / width_ratio)) + 1
    num_x = np.int32(n / num_y) + 1
    
    # create regularly spaced neurons
    x = np.linspace(0., shape[1]-1, num_x, dtype=np.float32)
    y = np.linspace(0., shape[0]-1, num_y, dtype=np.float32)
    coords = np.stack(np.meshgrid(x, y), -1).reshape(-1,2)
    
    # compute spacing
    init_dist = np.min((x[1]-x[0], y[1]-y[0]))
    min_dist = init_dist * (1 - sensitivity)
    
    assert init_dist >= min_dist
    print(min_dist)
    
    # perturb points
    max_movement = (init_dist - min_dist)/2
    noise = np.random.uniform(
        low=-max_movement,
        high=max_movement,
        size=(len(coords), 2))
    coords += noise
    
    # plot
    plt.figure(figsize=(10*width_ratio,10))
    plt.scatter(coords[:,0], coords[:,1], s=3)
    plt.show()
    

    【讨论】:

      【解决方案3】:

      据我了解,您正在寻找一种算法来在一个框中创建许多随机点,这样没有两个点比某个最小距离更近。如果这是您的问题,那么您可以利用统计物理学,并使用分子动力学软件解决它。此外,您确实需要分子动力学或蒙特卡罗来获得该问题的精确解。

      您将 N 个原子放在一个矩形框中,在它们之间创建一个固定半径的排斥相互作用(例如移位的 Lennard-Jones 相互作用),然后运行一段时间的模拟(直到您看到这些点在整个盒子)。根据统计物理定律,您可以证明点的位置将是最大随机,前提是点不能接近某个距离。 如果您使用迭代算法,例如将点一个接一个放置并在它们重叠时拒绝它们,这将是不正确的

      我估计 10000 点的运行时间为几秒钟,100k 的运行时间为几分钟。我使用 OpenMM 进行所有分子动力学模拟。

      【讨论】:

        【解决方案4】:
        #example of generating 50 points in a square of 4000x4000 and with minimum distance of 400
        
        import numpy as np
        import random as rnd
        
        n_points=50
        x,y = np.zeros(n_points),np.zeros(n_points)
        x[0],y[0]=np.round(rnd.uniform(0,4000)),np.round(rnd.uniform(0,4000))
        min_distances=[]
        i=1
        while i<n_points :
            x_temp,y_temp=np.round(rnd.uniform(0,4000)),np.round(rnd.uniform(0,4000))
            distances = []
            for j in range(0,i):
                distances.append(np.sqrt((x_temp-x[j])**2+(y_temp-y[j])**2))
            min_distance = np.min(distances)
            if min_distance>400 :
                min_distances.append(min_distance)
                x[i]=x_temp
                y[i]=y_temp
                i = i+1
        
        print(x,y)
        

        【讨论】:

          猜你喜欢
          • 2021-12-27
          • 1970-01-01
          • 2022-01-18
          • 2019-05-08
          • 1970-01-01
          • 1970-01-01
          • 1970-01-01
          • 2018-05-25
          • 2013-11-05
          相关资源
          最近更新 更多