Python/numpy棘手的切片问题 [英] Python/numpy tricky slicing problem
问题描述
我遇到一些麻木的问题.我需要一个numpy数组，通过返回一个切片作为我切片的数据的视图(而不是副本)来以一种不寻常的方式运行.所以这是我想做的一个例子:
说我们有一个像这样的简单数组:
a = array([1, 0, 0, 0])
我想使用以下语法，使用数组中的上一个条目更新数组中的连续条目(从左向右移动):
a[1:] = a[0:3]
这将得到以下结果:
a = array([1, 1, 1, 1])
或者类似这样的东西:
a[1:] = 2*a[:3]
# a = [1,2,4,8]
为进一步说明，我需要以下行为:
for i in range(len(a)):
if i == 0 or i+1 == len(a): continue
a[i+1] = a[i]
除了我想要numpy的速度.
numpy的默认行为是获取切片的副本，所以我实际得到的是:
a = array([1, 1, 0, 0])
我已经将此数组作为ndarray的子类，因此如果需要我可以对其进行进一步的更改，我只需要不断更新右侧的切片，因为它会更新左侧的切片边.
我在做梦吗?或者这可能是魔术吗?
更新:这都是因为我试图或多或少地使用Gauss-Seidel迭代来解决线性代数问题.这是一个涉及谐波函数的特例，我试图避免涉及这一点，因为它确实不是必需的，并且可能会使事情进一步混淆，但这是可行的.
算法是这样的:
while not converged:
for i in range(len(u[:,0])):
for j in range(len(u[0,:])):
# skip over boundary entries, i,j == 0 or len(u)
u[i,j] = 0.25*(u[i-1,j] + u[i+1,j] + u[i, j-1] + u[i,j+1])
对吗?但是，您可以通过两种方式进行操作，Jacobi涉及与相邻元素一起更新每个元素，而无需考虑在while循环循环之前已经进行的更新，要在循环中执行此操作，您需要复制数组，然后从复制的数组中更新一个数组.但是，Gauss-Seidel使用您已经为i-1和j-1条目中的每一个更新的信息，因此不需要复制，因为每次更新单个元素后都会对数组进行重新评估，所以循环本质上应该是知道"的.也就是说，每当我们调用u [i-1，j]或u [i，j-1]之类的条目时，在上一个循环中计算出的信息就会在那里.
我想用numpy切片用一段漂亮的干净代码替换这种缓慢而丑陋的嵌套循环情况:
u[1:-1,1:-1] = 0.25(u[:-2,1:-1] + u[2:,1:-1] + u[1:-1,:-2] + u[1:-1,2:])
但是结果是Jacobi迭代，因为当您获取一个切片:u [:，-2,1:-1]时，您将复制数据，因此该切片不知道进行了任何更新.现在numpy仍然循环吗?它不是并行的，只是一种更快的循环方式，看起来像python中的并行操作.我想通过窃取numpy来利用此行为，以便在获取切片时返回指针而不是副本.正确的?然后，每次numpy循环时，该片将更新"或实际上只是复制更新中发生的一切.为此，我需要在数组的两侧都将切片作为指针.
无论如何，如果那里真的有一个非常聪明的人真是太棒了，但是我已经辞职了，以为唯一的答案就是循环使用C.
最新答案，但这在Google上出现了，所以我可能指向OP想要的文档.您的问题很清楚:使用NumPy切片时，会创建临时对象.通过快速调用weave.blitz来包装代码，以摆脱临时工作，并具有您想要的行为.
有关完整的详细信息，请阅读 PerformancePython教程的weave.blitz部分.. >
I have a problem with some numpy stuff. I need a numpy array to behave in an unusual manner by returning a slice as a view of the data I have sliced, not a copy. So heres an example of what I want to do:
Say we have a simple array like this:
a = array([1, 0, 0, 0])
I would like to update consecutive entries in the array (moving left to right) with the previous entry from the array, using syntax like this:
a[1:] = a[0:3]
This would get the following result:
a = array([1, 1, 1, 1])
Or something like this:
a[1:] = 2*a[:3]
# a = [1,2,4,8]
To illustrate further I want the following kind of behaviour:
for i in range(len(a)):
if i == 0 or i+1 == len(a): continue
a[i+1] = a[i]
Except I want the speed of numpy.
The default behavior of numpy is to take a copy of the slice, so what I actually get is this:
a = array([1, 1, 0, 0])
I already have this array as a subclass of the ndarray, so I can make further changes to it if need be, I just need the slice on the right hand side to be continually updated as it updates the slice on the left hand side.
Am I dreaming or is this magic possible?
Update: This is all because I am trying to use Gauss-Seidel iteration to solve a linear algebra problem, more or less. It is a special case involving harmonic functions, I was trying to avoid going into this because its really not necessary and likely to confuse things further, but here goes.
The algorithm is this:
while not converged:
for i in range(len(u[:,0])):
for j in range(len(u[0,:])):
# skip over boundary entries, i,j == 0 or len(u)
u[i,j] = 0.25*(u[i-1,j] + u[i+1,j] + u[i, j-1] + u[i,j+1])
Right? But you can do this two ways, Jacobi involves updating each element with its neighbours without considering updates you have already made until the while loop cycles, to do it in loops you would copy the array then update one array from the copied array. However Gauss-Seidel uses information you have already updated for each of the i-1 and j-1 entries, thus no need for a copy, the loop should essentially 'know' since the array has been re-evaluated after each single element update. That is to say, every time we call up an entry like u[i-1,j] or u[i,j-1] the information calculated in the previous loop will be there.
I want to replace this slow and ugly nested loop situation with one nice clean line of code using numpy slicing:
u[1:-1,1:-1] = 0.25(u[:-2,1:-1] + u[2:,1:-1] + u[1:-1,:-2] + u[1:-1,2:])
But the result is Jacobi iteration because when you take a slice: u[:,-2,1:-1] you copy the data, thus the slice is not aware of any updates made. Now numpy still loops right? Its not parallel its just a faster way to loop that looks like a parallel operation in python. I want to exploit this behaviour by sort of hacking numpy to return a pointer instead of a copy when I take a slice. Right? Then every time numpy loops, that slice will 'update' or really just replicate whatever happened in the update. To do this I need slices on both sides of the array to be pointers.
Anyway if there is some really really clever person out there that awesome, but I've pretty much resigned myself to believing the only answer is to loop in C.
Late answer, but this turned up on Google so I probably point to the doc the OP wanted. Your problem is clear: when using NumPy slices, temporaries are created. Wrap your code in a quick call to weave.blitz to get rid of the temporaries and have the behaviour your want.
Read the weave.blitz section of PerformancePython tutorial for full details.
这篇关于Python/numpy棘手的切片问题的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！
