On 17/10/06, Travis Oliphant <[EMAIL PROTECTED]> wrote:
Charles R Harris wrote:
>
> In [3]: array([[1,2,3,4,5,6]]).reshape((3,2),order='F')
> Out[3]:
> array([[1, 4],
> [2, 5],
> [3, 6]])
>
> I also don't understand why a copy is returned if 'F' just fiddles
> with the indices and strides; the underlying data should be the same,
> just the view changes. FWIW, I think both examples should be returning
> views.
You are right, it doesn't need to. My check is not general enough.
It can be challenging to come up with a general way to differentiate the
view-vs-copy situation and I struggled with it. In this case, it's the
fact that while self->nd > 1, the other dimensions are only of shape 1
and so don't really matter. If you could come up with some kind of
striding check that would distinguish the two cases, I would appreciate it.
I wrote, just for exercise I guess, a routine that checks to see if a
copy is needed for a reshape (and if not, it computes the resulting
strides) for an arbitrary array. (It currently only does a C-order
reshape, but F-order would be an easy addition.) It does, however,
successfully handle arbitrarily strided arrays with arbitrary
arrangements of "1" dimensions having arbitrary strides. I think it
also correctly handles 0 strides. I don't imagine you'd want to *use*
it, but it does provide a complete solution to copy-less reshaping.
I'd be happy to help with the C implementation built into numpy, but
for the life of me I can't figure out how it works.
A. M. Archibald
P.S. I haven't written down a *proof* that this implementation never
copies unnecessarily, but I'm pretty sure that it doesn't. It has a
reshape_restricted that effectively does a ravel() first; this can't
cope with things like length (5,2,3) -> length (5,3,2) unless the
strides are conveniently equal, so there's a wrapper, reshape(), that
breaks the reshape operation up into pieces that can be done by
reshape_restricted. -AMA
from numpy import multiply, array, empty, zeros, add, asarray
from itertools import izip
class MustCopy(Exception):
pass
def index(indices, strides):
indices = asarray(indices)
strides = asarray(strides)
if len(indices)!=len(strides):
raise ValueError
return add.reduce(indices*strides)
def reshape_restricted(from_strides,from_lengths,to_lengths):
if multiply.reduce(from_lengths)!=multiply.reduce(to_lengths):
raise ValueError, "Array sizes are not compatible"
for i in xrange(1,len(from_strides)):
if not (from_strides[i-1]==from_strides[i]*from_lengths[i] or
from_lengths[i]==1): # If a dimension has length one, who cares what its stride is?
raise MustCopy
r = empty(len(to_lengths), dtype=int)
j = len(from_lengths)-1
while j>=0 and from_lengths[j]==1:
j-=1 # Ignore strides on length-0 dimensions
if j<0: # Both arrays have size 1
r[-1] = 1
else:
r[-1] = from_strides[j]
i = -1
while len(to_lengths)+i>0:
r[i-1]=to_lengths[i]*r[i]
i -= 1
return r
def reshape(from_strides,from_lengths,to_lengths):
if multiply.reduce(from_lengths)!=multiply.reduce(to_lengths):
raise ValueError
r = empty(len(to_lengths), dtype=int)
j1 = len(from_lengths)
j2 = len(to_lengths)
while j1!=0:
i1 = j1-1
i2 = j2-1
p1 = from_lengths[i1]
p2 = to_lengths[i2]
while p1!=p2:
if p1<p2:
i1-=1
assert i1>=0
p1*=from_lengths[i1]
else:
i2-=1
assert i2>=0
p2*=to_lengths[i2]
r[i2:j2] = reshape_restricted(from_strides[i1:j1],
from_lengths[i1:j1],
to_lengths[i2:j2])
j1 = i1
j2 = i2
return r
def test_reshape(from_strides,from_lengths,to_strides,to_lengths):
if multiply.reduce(from_lengths)!=multiply.reduce(to_lengths):
raise ValueError
if len(from_strides)!=len(from_lengths) or len(to_strides)!=len(to_lengths):
raise ValueError
def walk(lengths):
i = zeros(len(lengths),int)
while True:
yield i
j = len(lengths)-1
while i[j]==lengths[j]-1:
i[j]=0
j-=1
if j<0:
return
i[j] += 1
for (i1,i2) in izip(walk(from_lengths), walk(to_lengths)):
if index(i1,from_strides)!=index(i2,to_strides):
raise ValueError
-------------------------------------------------------------------------
Using Tomcat but need to do more? Need to support web services, security?
Get stuff done quickly with pre-integrated technology to make your job easier
Download IBM WebSphere Application Server v.1.0.1 based on Apache Geronimo
http://sel.as-us.falkag.net/sel?cmd=lnk&kid=120709&bid=263057&dat=121642
_______________________________________________
Numpy-discussion mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/numpy-discussion