@ Теория Пола совершенно верна.В этом ответе я использую perf и отладчик для погружения, чтобы поддержать эту теорию.
Во-первых, давайте посмотрим, где тратится время выполнения (точные сведения см. В листингах для run.py ниже).код).
Для n=1 мы видим следующее:
Event count (approx.): 3388750000
Overhead Command Shared Object Symbol
34,04% python umath.cpython-36m-x86_64-linux-gnu.so [.] DOUBLE_less
32,71% python multiarray.cpython-36m-x86_64-linux-gnu.so [.] _aligned_strided_to_contig_size8_srcstride0
28,16% python libc-2.23.so [.] __memmove_ssse3_back
1,46% python multiarray.cpython-36m-x86_64-linux-gnu.so [.] PyArray_TransferNDimToStrided
по сравнению с n=2:
Event count (approx.): 28954250000
Overhead Command Shared Object Symbol
40,85% python libc-2.23.so [.] __memmove_ssse3_back
40,16% python multiarray.cpython-36m-x86_64-linux-gnu.so [.] PyArray_TransferNDimToStrided
8,61% python umath.cpython-36m-x86_64-linux-gnu.so [.] DOUBLE_less
8,41% python multiarray.cpython-36m-x86_64-linux-gnu.so [.] _contig_to_contig
Для n = 2 событий в 8,5 раз большесчитается, но только для удвоенных данных, поэтому нам нужно объяснить фактор замедления 4.
Еще одно важное наблюдение: во время выполнения преобладают операции с памятью для n=2 и (менее очевидно) такжедля n=1 (_aligned_strided_to_contig_size8_srcstride0 это все о копировании данных), они завышают затраты на сравнение - DOUBLE_less.
Очевидно, что PyArray_TransferNDimtoStrided вызывается для обоих размеров, поэтомупочему так велика разница в его доле времени выполнения?
Показанное собственное время PyArray_TransferNDimtoStrided - это не время, необходимое для копирования, а накладные расходы: указатели настраиваются так, чтобыв последнем измерении может быть скопирован за один раз с помощью stransfer:
PyArray_TransferNDimToStrided(npy_intp ndim,
....
/* A loop for dimensions 0 and 1 */
for (i = 0; i < shape1; ++i) {
if (shape0 >= count) {
stransfer(dst, dst_stride, src, src_stride0,
count, src_itemsize, data);
return 0;
}
else {
stransfer(dst, dst_stride, src, src_stride0,
shape0, src_itemsize, data);
}
count -= shape0;
src += src_stride1;
dst += shape0*dst_stride;
}
...
Эти функции передачи являются _aligned_strided_to_contig_size8_srcstride0 (см. сгенерированный код вперечисление ниже) и _contig_to_contig:
_contig_to_contig используется в случае n=2 и tranfers 2-double (последнее измерение имеет 2 значения), накладные расходы на настройку указателей довольно высоки! _aligned_strided_to_contig_size8_srcstride0 используется для n=1 и передает 1000 дубликатов за вызов (как указал @Paul икак мы скоро увидим, numpy достаточно умен, чтобы отбросить размеры, которые состоят из 1 элемента), накладными расходами на настройку указателей можно пренебречь.
Кстати, эти функции используются вместопростой цикл for для использования векторизации современных процессоров: с шагом, известным во время компиляции, компилятор может векторизовать код (что компиляторы часто не могут сделать для шагов, известных только во время выполнения), таким образом, numpy анализирует схему доступаи отправляет в разные предварительно скомпилированные функции.
Остался один вопрос: действительно ли numpy отбрасывает последнее измерение, если его размер равен 1, как показывают наши наблюдения?
Это легко проверить с помощью отладчика:
Asдля коэффициента скорости 4, который «теряется» при сравнении n=2 с n=1: он не имеет особого значения и является лишь случайным значением на моей машине: изменение размера матрицы с 10 ^ 3 до 10^ 4 еще больше сместит преимущество (меньше накладных расходов) на n=1, что приведет к тому, что на моей машине будет коэффициент потери скорости 12.
run.py
import sys
import numpy as np
n=int(sys.argv[1])
x, y = (np.random.uniform(size=(1, 1000, n)),
np.random.uniform(size=(1000, 1, n)))
for _ in range(10000):
y<x
, а затем:
perf record python run.py 1
perf report
....
perf record python run.py 2
perf report
Сгенерированный источник _aligned_strided_to_contig_size8_srcstride0:
/*
* specialized copy and swap for source stride 0,
* interestingly unrolling here is like above is only marginally profitable for
* small types and detrimental for >= 8byte moves on x86
* but it profits from vectorization enabled with -O3
*/
#if (0 == 0) && 1
static NPY_GCC_OPT_3 void
_aligned_strided_to_contig_size8_srcstride0(char *dst,
npy_intp dst_stride,
char *src, npy_intp NPY_UNUSED(src_stride),
npy_intp N, npy_intp NPY_UNUSED(src_itemsize),
NpyAuxData *NPY_UNUSED(data))
{
#if 8 != 16
# if !(8 == 1 && 1)
npy_uint64 temp;
# endif
#else
npy_uint64 temp0, temp1;
#endif
if (N == 0) {
return;
}
#if 1 && 8 != 16
/* sanity check */
assert(npy_is_aligned(dst, _ALIGN(npy_uint64)));
assert(npy_is_aligned(src, _ALIGN(npy_uint64)));
#endif
#if 8 == 1 && 1
memset(dst, *src, N);
#else
# if 8 != 16
temp = _NPY_NOP8(*((npy_uint64 *)src));
# else
# if 0 == 0
temp0 = (*((npy_uint64 *)src));
temp1 = (*((npy_uint64 *)src + 1));
# elif 0 == 1
temp0 = _NPY_SWAP8(*((npy_uint64 *)src + 1));
temp1 = _NPY_SWAP8(*((npy_uint64 *)src));
# elif 0 == 2
temp0 = _NPY_SWAP8(*((npy_uint64 *)src));
temp1 = _NPY_SWAP8(*((npy_uint64 *)src + 1));
# endif
# endif
while (N > 0) {
# if 8 != 16
*((npy_uint64 *)dst) = temp;
# else
*((npy_uint64 *)dst) = temp0;
*((npy_uint64 *)dst + 1) = temp1;
# endif
# if 1
dst += 8;
# else
dst += dst_stride;
# endif
--N;
}
#endif/* @elsize == 1 && 1 -- else */
}
#endif/* (0 == 0) && 1 */