时间戳

首页 > 解决方案 > 如何使用 numpy 执行此复制操作？

python - 如何使用 numpy 执行此复制操作？

问题描述

我一直在研究“扩散蒙特卡罗”的基本模拟，以找到氢分子的基态能量。算法的一个关键部分正在痛苦地减慢我的代码，我不知道如何修复它。

这就是代码正在做的事情。我有一个名为 x 的 6 x N numpy 数组。该阵列代表 N 个随机“步行者”，它们对 6 维相空间进行采样（两个电子乘以 3 维是 6 维）。我建议对每个“walker”进行某些随机更改以获得我的新“walker”，然后使用公式为每个新 walker 吐出一个数字“m”。

数字 m 可以是 0、1、2 或 3。这是困难的部分。如果 m 为 0，则它对应的“walker”将从数组中删除。如果 m 为 1，则步行者仍然存在。如果 m 是 2，那么 walker 仍然存在，我必须在数组中创建 walker 的新副本。如果 m 为 3，则步行者仍然存在，我必须在数组中制作步行者的两个新副本。此后代码重复；对我的步行者数组等提出了随机更改。

所以; 以下是减慢算法的代码。这是最后一部分的代码，我必须通过我的 m 并确定如何处理每个“walker”，并创建我的新数组 x 以用于算法的下一次迭代。

        j1 = 0
        n1 = len(x[0,:])
        x_n = np.ones((6,2*n1))
        for i in range(0,n1):
            if m[i] == 1:
                x_n[:,j1] = x[:,i]
                j1 = j1 + 1
            if m[i] == 2:
                x_n[:,j1] = x[:,i]
                x_n[:,j1+1] = x[:,i]
                j1 = j1 + 2
            if m[i] == 3:
                x_n[:,j1] = x[:,i]
                x_n[:,j1+1] = x[:,i]
                j1 = j1 + 3
        x = np.ones((6,j1))
        for j in range(0,j1):
            x[:,j] = x_n[:,j]

我的问题如下；有没有办法使用 numpy 本身来做我在这段代码中所做的事情？根据我的经验，Numpy 往往比 for 循环快得多。在变分蒙特卡罗模拟中直接使用 numpy，我能够在运行时间上实现 100 倍的改进。如果您想要完整的代码来实际运行算法，那么我可以发布它；它只是相当长。

标签： pythonnumpyphysics

解决方案

令 M 为每个随机游走者的 m 个值的 N x 1 数组。

让 X 成为您原来的 6 x N 数据数组

# np.where returns a list of indices where the condition is satisfied 
zeros =  np.where(M == 0) # don't actually need this variable, I just did it for completeness
ones =   np.where(M == 1)
twos =   np.where(M == 2)
threes = np.where(M == 3)

# use the lists of indices to access the relevant portions of X
ones_array = X[:,ones]
twos_array = X[:,twos]
threes_array = X[:,threes]

# update X with one copy where m = 1, two copies where m = 2, three copies where m = 3
X = np.concatenate((ones_array,twos_array,twos_array,threes_array,threes_array,threes_array),axis = 1)

这不会保留步行者的顺序，所以如果这很重要，代码会略有不同。

推荐阅读