[Numpy-discussion] storage order and concatenate

[David Cournapeau]

Travis Oliphant wrote:
> David Cournapeau wrote:
>   
>> Hi,
>>
>>
>> What are the rules concerning storage with numpy ? 
>>     
> The rule is that a numpy array has "strides" which specify how many 
> bytes to skip to get to the next element in the array.   That's the 
> internal model.  There are no hard and fast rules about storage order.  
> Internally, C-order is as good as Fortran-order (except the array 
> iterator gives special preference to C-order and all functions for which 
> the order can be specified (like zeros) default to C-order).
>   
Ok, this is again a bad habit from matlab to think in C or F order...
> Thus, the storage order is whatever the strides say it is.  Now, there 
> are flags that keep track of whether or not the strides agree with the 2 
> recognized special cases of "Fortran-order" (first-index varies the 
> fastest) or "C-order" (last-index varies the fastest).  But, this is 
> only for convenience.   Very few functions actually require a 
> specification of storage order.  Those that allow it default to "C-order".
>
> You can't think of a NumPy array has having a particular storage order 
> unless you explicitly request it.  One of the most common ways that 
> Fortran-order arrays show up, for example is when a C-order array is 
> transposed.  A transpose operation does nothing except flip the strides 
> (and therefore the flags) of the array.  This is what is happening in 
> concatenate (using axis=1) to give you a Fortran-order array.  
> Bascially, code equivalent to the following is being run:  
> concatenate([X1.T, X2.T]).T
>
> In the second example, you explicitly create the array (and therefore 
> the strides) as C-order and then fill it (so it doesn't change on you).  
> The first example used array calculations which don't guarantee the 
> storage order. 
>
> This is all seamless to the user until you have to interface with 
> extension code.  Ideally, you write compiled code that deals with 
> strided arrays.  If you can't, then you request an array of the required 
> storage-order. 
>
> By the way, for interfacing with ctypes, check out the 
> ctypeslib.ndpointer class-creation function for flag checking and the 
> require function for automatic conversion of an array to specific 
> requirements.
>
>   
I tried to to that at first, but I couldn't make the examples of numpy 
works; after having tried at home with beta versions, it looks like it 
is a ctype version problem, as it works with ctype 1.0 + numpy 1.0rc1.

Thanks for the explanations, this answers most questions I had in mind 
for numpy internal layout compared to matlab which I am used to. I think 
this should be in the wiki somewhere; do you mind if I use your email as 
a basis for the tentative numpy tutorial (memory layout section, maybe) ?

David