Addition speed of numpy arrays with different contiguous-type
<p>Numpy arrays are stored with different contiguous types (C- and F-). When using numpy.swapaxes(), the contiguous type gets changed. I need to add two multidimensional arrays (3d to be more specific), one of which comes from another array with swapped axes. What I've noticed is that when the first axis gets swapped with the last axis, in the case of a 3d array, the contiguous type changes from C- to F-. And adding two arrays with different contiguous type is extremely slow (~6 times slower than adding two C-contiguous arrays). However, if other axes are swapped (0-1 or 1-2), the resulting array would have false flags for both C- and F- contiguous (non-contiguous). The weird thing to me is that adding one array of C-configuous and one array neither C- nor F- contiguous, is in fact only slightly slower than adding two arrays of same type. Here are my two questions:</p>

<ol>
<li><p>Why does it seem to be different for C-&amp;F-contiguous arrray addition and C-&amp;non-contiguous array addition? Is is caused by different rearranging mechanism or simply because the rearranging distance between C- and F- contiguous is longest for all possible axes orders?</p></li>
<li><p>If I have to add a C-contiguous array and a F-contiguous/non-contiguous array, what is the best way to accelerate the speed?</p></li>
</ol>

<p>Below is a minimum example of what I encountered. The three printed durations on my computer are 2.0s (C-contiguous + C-contiguous), 12.4s (C-contiguous + F-contiguous), 3.4s (C-contiguous + non-contiguous) and 3.3s (C-contiguous + non-contiguous).</p>

<pre><code>import numpy as np
import time

np.random.seed(1234)

a = np.random.random((300, 400, 500))  # C-contiguous
b = np.swapaxes(np.random.random((500, 400, 300)), 0, 2)  # F-contiguous
c = np.swapaxes(np.random.random((300, 500, 400)), 1, 2)  # Non-contiguous
d = np.swapaxes(np.random.random((400, 300, 500)), 0, 1)  # Non-contiguous

t = time.time()
for n in range(10):
    result = a + a
print(time.time() - t)

t = time.time()
for n in range(10):
    result = a + b
print(time.time() - t)

t = time.time()
for n in range(10):
    result = a + c
print(time.time() - t)

t = time.time()
for n in range(10):
    result = a + d
print(time.time() - t)
</code></pre>

<p>These types (<code>F</code> and <code>C</code>) denote whether a matrix (or multi-dimensional array) is stored in column-major (<code>C</code> as in C language which uses column-major storage) or row-major (<code>F</code> as in Fortran language which uses row-major storage).</p>

<p>Both do not really vary in speed. It is just a abstraction layer. No matter which one you use, it brings performance wise the same. </p>

<p>However, what makes an enormous difference is whether arrays are contiguous or not. If they are contiguous you will have good timings cause of caching effects, vectorization and other optimization games that the compiler might apply.</p>

python arrays performance numpy contiguous