Я не уверен, что кэширование действительно является единственным наиболее влиятельным фактором.
Я также не являюсь специалистом в области информатики, поэтому я вполне могу ошибаться, но позвольте мне провести вас через пару наблюдений. Для простоты я использую вызов @hpaulj, который «+» показывает практически тот же эффект, что и «**».
Моя рабочая гипотеза заключается в том, что это накладные расходы на внешние циклы, которые, как я считаю, существенно дороже, чем непрерывные векторизуемые самые внутренние циклы.
Итак, давайте сначала минимизируем количество повторяющихся данных, чтобы кэширование вряд ли оказало большое влияние:
>>> from timeit import repeat
>>> import numpy as np
>>>
>>> def mock_data(k, N, M):
... x = list(np.random.randint(0, 10000, (k, N, M)))
... y = list(np.random.randint(0, 10000, (k, M)))
... z = list(np.random.randint(0, 10000, (k, N, 1)))
... return x, y, z
...
>>> k, N, M = 500, 5000, 4
>>>
>>> repeat('x.pop() + y.pop()', setup='x, y, z = mock_data(k, M, N)', globals=globals(), number=k)
[0.017986663966439664, 0.018148145987652242, 0.018077059998176992]
>>> repeat('x.pop() + y.pop()', setup='x, y, z = mock_data(k, N, M)', globals=globals(), number=k)
[0.026680009090341628, 0.026304758968763053, 0.02680662798229605]
Здесь оба сценария имеют непрерывные данные, одинаковое количество дополнений, но версия с 5000 внешними итерациями существенно медленнее. Когда мы возвращаем кеширование, хотя и в нескольких испытаниях, разница остается примерно такой же, но соотношение становится еще более выраженным:
>>> repeat('x[0] + y[0]', setup='x, y, z = mock_data(k, M, N)', globals=globals(), number=k)
[0.011324503924697638, 0.011121788993477821, 0.01106808998156339]
>>> repeat('x[0] + y[0]', setup='x, y, z = mock_data(k, N, M)', globals=globals(), number=k)
[0.020170683041214943, 0.0202067659702152, 0.020624138065613806]
Возвращаясь к исходному сценарию «внешней суммы», мы видим, что в случае несмежных длинных измерений дело обстоит еще хуже. Поскольку нам нужно прочитать не больше данных, чем в непрерывном сценарии, это нельзя объяснить тем, что данные не кэшируются.
>>> repeat('z.pop() + y.pop()', setup='x, y, z = mock_data(k, M, N)', globals=globals(), number=k)
[0.013918839977122843, 0.01390116906259209, 0.013737019035033882]
>>> repeat('z.pop() + y.pop()', setup='x, y, z = mock_data(k, N, M)', globals=globals(), number=k)
[0.0335254140663892, 0.03351909795310348, 0.0335453050211072]
Кроме того, оба получают прибыль от пробного кэширования:
>>> repeat('z[0] + y[0]', setup='x, y, z = mock_data(k, M, N)', globals=globals(), number=k)
[0.012061356916092336, 0.012182610924355686, 0.012071475037373602]
>>> repeat('z[0] + y[0]', setup='x, y, z = mock_data(k, N, M)', globals=globals(), number=k)
[0.03265167598146945, 0.03277428599540144, 0.03247103898320347]
С точки зрения кашиста это в лучшем случае неубедительно.
Итак, давайте посмотрим на источник. После сборки текущего NumPy из архива вы найдете где-то в дереве почти 15000 строк сгенерированного компьютером кода в файле с именем «loops.c». Эти циклы являются самыми внутренними циклами ufuncs, наиболее релевантным для нашей ситуации является
#define BINARY_LOOP\
char *ip1 = args[0], *ip2 = args[1], *op1 = args[2];\
npy_intp is1 = steps[0], is2 = steps[1], os1 = steps[2];\
npy_intp n = dimensions[0];\
npy_intp i;\
for(i = 0; i < n; i++, ip1 += is1, ip2 += is2, op1 += os1)
/*
* loop with contiguous specialization
* op should be the code working on `tin in1`, `tin in2` and
* storing the result in `tout * out`
* combine with NPY_GCC_OPT_3 to allow autovectorization
* should only be used where its worthwhile to avoid code bloat
*/
#define BASE_BINARY_LOOP(tin, tout, op) \
BINARY_LOOP { \
const tin in1 = *(tin *)ip1; \
const tin in2 = *(tin *)ip2; \
tout * out = (tout *)op1; \
op; \
}
etc.
Полезная нагрузка в нашем случае кажется достаточно скудной, особенно если я правильно интерпретирую комментарий о непрерывной специализации и автовекторизации. Теперь, если мы делаем только 4 итерации, отношение накладных расходов к полезной нагрузке начинает выглядеть немного тревожно, и на этом оно не заканчивается.
В файле ufunc_object.c находим следующий фрагмент
/*
* If no trivial loop matched, an iterator is required to
* resolve broadcasting, etc
*/
NPY_UF_DBG_PRINT("iterator loop\n");
if (iterator_loop(ufunc, op, dtypes, order,
buffersize, arr_prep, arr_prep_args,
innerloop, innerloopdata) < 0) {
return -1;
}
return 0;
реальный цикл выглядит так
NPY_BEGIN_THREADS_NDITER(iter);
/* Execute the loop */
do {
NPY_UF_DBG_PRINT1("iterator loop count %d\n", (int)*count_ptr);
innerloop(dataptr, count_ptr, stride, innerloopdata);
} while (iternext(iter));
NPY_END_THREADS;
innerloop — это внутренний цикл, который мы рассмотрели выше. Сколько накладных расходов связано с iternext?
Для этого нам нужно обратиться к файлу nditer_templ.c.src, где мы находим
/*NUMPY_API
* Compute the specialized iteration function for an iterator
*
* If errmsg is non-NULL, it should point to a variable which will
* receive the error message, and no Python exception will be set.
* This is so that the function can be called from code not holding
* the GIL.
*/
NPY_NO_EXPORT NpyIter_IterNextFunc *
NpyIter_GetIterNext(NpyIter *iter, char **errmsg)
{
etc.
Эта функция возвращает указатель функции на одну из вещей, которую делает предварительная обработка.
/* Specialized iternext (@const_itflags@,@tag_ndim@,@tag_nop@) */
static int
npyiter_iternext_itflags@tag_itflags@_dims@tag_ndim@_iters@tag_nop@(
NpyIter *iter)
{
etc.
Разобрать это выше моих сил, но в любом случае это указатель на функцию, который должен вызываться на каждой итерации внешнего цикла, и, насколько мне известно, указатели на функции не могут быть встроены, поэтому по сравнению с 4 итерациями тривиального тела цикла это будет существенным.
Я, вероятно, должен профилировать это, но моих навыков недостаточно.
person
Paul Panzer
schedule
15.01.2018
+. - person hpaulj schedule 14.01.2018xиyявляются числами с плавающей запятой, относительные тайминги переключаются (m2быстрее, чемm1). - person Daniel F schedule 15.01.2018