【问题标题】:Rolling/Increasing dimensionality of a NumPy arrayNumPy 数组的滚动/增加维度
【发布时间】:2019-01-06 07:27:48
【问题描述】:

我目前正在尝试找到一种简单的方法来对 Python 中的 N 维数组执行以下操作。为简单起见,让我们从大小为 4 的一维数组开始。

X = np.array([1,2,3,4])

我要做的是创建一个新数组,将其命名为 Y,这样:

Y = np.array([1,2,3,4],[2,3,4,1],[3,4,1,2],[4,1,2,3])

所以我要做的是创建一个数组 Y,这样:

Y[:,i] = np.roll(X[:],-i, axis = 0)

我知道如何使用 for 循环来做到这一点,但我正在寻找一种更快的方法。我试图这样做的实际数组是一个 3 维数组,称之为 X。我正在寻找的是一种查找数组 Y 的方法,例如:

Y[:,:,:,i,j,k] = np.roll(X[:,:,:],(-i,-j,-k),axis = (0,1,2))

我可以通过使用 for 循环的 itertools.product 类来做到这一点,但这很慢。如果有人有更好的方法,请告诉我。我还用 GTX-970 安装了 CUPY,所以如果有办法使用 CUDA 更快地做到这一点,请告诉我。如果有人想要更多上下文,请告诉我。

这是我计算位置空间两点相关函数的原始代码。数组 x0 是一个 n × n × n 实数值数组,表示一个实标量字段。函数 iterate(j,s) 运行 j 次迭代。每次迭代都包括为每个格点生成 -s 和 s 之间的随机浮点数。然后计算动作 dS 的变化,并以 min(1,exp^(-dS)) 的概率接受变化

def momentum(k,j,s):
global Gxa
Gx = numpy.zeros((n,n,t))
for i1 in range(0,k):
    iterate(j,s)
    for i2,i3,i4 in itertools.product(range(0,n),range(0,n),range(0,n)):
        x1 = numpy.roll(numpy.roll(numpy.roll(x0, -i2, axis = 0),-i3, axis = 1),-i4,axis = 2)
        x2 = numpy.mean(numpy.multiply(x0,x1))
        Gx[i2,i3,i4] = x2
    Gxa = Gxa + Gx
Gxa = Gxa/k

【问题讨论】:

    标签: python arrays numpy numpy-ndarray


    【解决方案1】:

    方法#1

    我们可以在这里将this idea 扩展到我们的3D 数组案例。因此,只需沿三个暗淡连接切片版本,然后使用基于scikit-image's view_as_windowsnp.lib.stride_tricks.as_strided 有效地获得最终输出作为连接版本的跨步视图,就像这样 -

    from skimage.util.shape import view_as_windows
    
    X1 = np.concatenate((X,X[:,:,:-1]),axis=2)
    X2 = np.concatenate((X1,X1[:,:-1,:]),axis=1)
    X3 = np.concatenate((X2,X2[:-1,:,:]),axis=0)
    out = view_as_windows(X3,X.shape)
    

    方法#2

    对于非常大的数组,我们可能想要初始化输出数组,然后重新使用之前方法中的X3 来分配它并对其进行切片。这种切片过程将比原始滚动更快。实施将是 -

    m,n,r = X.shape
    Yout = np.empty((m,n,r,m,n,r),dtype=X.dtype)
    for i in range(m):
        for j in range(n):
            for k in range(r):
                Yout[:,:,:,i,j,k] = X3[i:i+m,j:j+n,k:k+r]
    

    【讨论】:

    • 太棒了!我得到了很大的加速。坏消息;我现在在三个维度上限制为大约 33 的最大晶格大小(我正在模拟晶格量子场论),因为如果超过 33^6 会出现内存错误,因为 33^6 是一个非常大的数字。有什么建议的方法来解决这个问题吗?
    • 它真的不起作用。我认为实际的内存接收器将是当我正在执行 np.multiply 以元素方式将新的六维矩阵与原始矩阵相乘时,因为我正在尝试为交互理论计算欧几里德(虚时间)绿色函数格子。无论使用您的方法,我都可以探测以前对我来说不可行的晶格尺寸,所以谢谢!
    • @chuckstables By - It doesn't really work. ,你是说方法 2 不起作用吗?如果是这样,您的意思是它的速度慢还是输出错误?我已经检查了性能和正确性。
    • 我的意思是方法 2。它给出了正确的输出,并且似乎与第一种方法一样快,但我仍然遇到内存限制。有时间我会再做一些测试,然后告诉你我发现了什么。
    • @chuckstables 只是好奇,您的滚动方法是否也达到了内存限制?如果是这样,我想我们没有其他办法了。
    猜你喜欢
    • 2012-11-12
    • 2021-03-20
    • 2019-12-19
    • 2015-05-26
    • 2011-10-12
    • 2011-03-04
    • 2015-01-29
    • 2017-06-10
    相关资源
    最近更新 更多