Есть ли более быстрый способ найти максимальное/минимальное значение в массиве, чем в OS X vecLib?

Question

Мне нужно найти самое высокое и самое низкое значение в массиве поплавков. Необязательно, я хочу, чтобы иметь возможность пропускать элементы массива и оценивать только каждую 2-ю, 4-й, 8-й и т.д. элемент:

float maxValue = 0.0; 
float minValue = 0.0; 

int index = 0; 
int stepwith = 8; 

for(int i = 0; i < pixelCount; i++) 
{ 
    float value = data[index]; 

    if(value > maxValue) 
      maxValue = value; 

    if(value < minValue) 
      minValue = value; 

    index += stepwidth; 
    if(index >= dataLength) 
     break; 
} 

Это, кажется, самый быстрый способ без использования другой магии.

Так что я попробовал другую магию, а именно vIsmax() и vIsmin() функцию из vecLib (включено в OSX»Ускорить рамки) которым, по-видимому использует ускорение процессора уровня векторных операций:

int maxIndex = vIsmax(pixelCount, data); 
int minIndex = vIsmin(pixelCount, data); 

float maxValue = data[maxIndex]; 
float minValue = data[minIndex]; 

Это очень быстро, но не позволяет пропускать значения (функции не предлагают аргумент «шаг»). Это делает его на самом деле медленнее, чем мой первый код, потому что я могу легко пропустить каждое восьмое значение.

Я даже нашел третий путь с BLAS которым реализует подобную функцию:

cblas_isamax(const int __N, const float *__X, const int __incX) 

с __incX = походкой перескочить значения, но это не так быстро, на всех, и только находит абсолютные максимумы которым Безразлично Я работаю для меня.

Итак, мой вопрос: может ли кто-нибудь подумать о другом способе ускорить это?

Answer 1

Следуя предложению в комментарии, я просмотрел Intel intrinsics и придумал этот код. Справедливое предупреждение: Это мой самый первый подход к внутренним характеристикам и очень экспериментальный. Он работает, хотя:

#include <emmintrin.h> 

void vec_minmax(float * inData, int length, float * outMin, float * outMax) 
{ 
    // In each iteration of the loop we will gather 8 from 64 
    // values (only fetching every 8th value). 

    // Build an index set that points to 8 consecutive floats. 
    // These indexes will later be spread up by factor 8 so they 
    // point to every 8th float. 
    // NOTE: these indexes are bytes, in reverse order. 
    __m256i vindex = _mm256_set_epi32(28, 24, 20, 16, 12, 8, 4, 0); 

    // Gather the first 8 floats. 
    __m256 v_min = _mm256_i32gather_ps(inData, vindex, 8); 
    __m256 v_max = v_min; 

    for(int i = 64; i < length; i += 64) 
    { 
     // gather the next set of floats. 
     __m256 v_cur = _mm256_i32gather_ps((inData + i), vindex, 8); 

     // Compare every member and store the results in v_min and v_max. 
     v_min = _mm256_min_ps(v_min, v_cur); 
     v_max = _mm256_max_ps(v_max, v_cur); 
    } 

    // Store the final result in two arrays. 
    float * max8; 
    float * min8; 

    posix_memalign((void **)&min8, 32, (8 * sizeof(float))); 
    posix_memalign((void **)&max8, 32, (8 * sizeof(float))); 

    _mm256_store_ps(min8, v_min); 
    _mm256_store_ps(max8, v_max); 

    // Find the min/max value in the arrays. 
    * outMin = min8[0]; 
    * outMax = max8[0]; 
    for(int i = 0; i < 8; i++) 
    { 
     if(min8[i] < * outMin) 
      * outMin = min8[i]; 

     if(max8[i] > * outMax) 
      * outMax = max8[i]; 
    } 
} 

Так эта функция находит минимальное и максимальное значение в наборе поплавков, рассматривая только каждый 8-значение, которое достаточно точности для моих потребностей.

К сожалению, это не значительно быстрее, чем тривиальный скалярный подход с простым for-loop и двумя if-операторами (как показано выше). По крайней мере, не с шагом 8.

Answer 2

Вот реализация для случая stride = 8 с использованием SSE. Я протестировал код, но у меня еще не было времени проверить его. Я не совсем уверен, что это будет быстрее, чем скалярная реализация, но это стоит попробовать ...

#include <tmmintrin.h> 
#include <float.h> 

void vec_minmax_stride_8(const float * inData, int length, float * outMin, float * outMax) 
{ 
    __m128i vmax = _mm_set1_ps(-FLT_MAX); 
    __m128i vmin = _mm_set1_ps(FLT_MAX); 

    for (int i = 0; i < length; i += 32) 
    { 
     __m128i v0 = _mm_loadu_ps(&inData[i]); 
     __m128i v1 = _mm_loadu_ps(&inData[i + 8]); 
     __m128i v2 = _mm_loadu_ps(&inData[i + 16]); 
     __m128i v3 = _mm_loadu_ps(&inData[i + 24]); 

     v0 = _mm_shuffle_ps(v0, v1, _MM_SHUFFLE(0, 0, 0, 0)); 
     v2 = _mm_shuffle_ps(v2, v3, _MM_SHUFFLE(0, 0, 0, 0)); 
     v0 = _mm_shuffle_ps(v0, v2, _MM_SHUFFLE(2, 0, 2, 0)); 

     vmax = _mm_max_ps(vmax, v0); 
     vmin = _mm_min_ps(vmin, v0); 
    } 

    vmax = _mm_max_ps(vmax, _mm_alignr_epi8(vmax, vmax, 4)); 
    vmin = _mm_min_ps(vmin, _mm_alignr_epi8(vmin, vmin, 4)); 

    vmax = _mm_max_ps(vmax, _mm_alignr_epi8(vmax, vmax, 8)); 
    vmin = _mm_min_ps(vmin, _mm_alignr_epi8(vmin, vmin, 8)); 

    _mm_store_ss(outMax, vmax); 
    _mm_store_ss(outMin, vmin); 
} 

Answer 3

исключения случаев, когда объем вычислений высок - этот пример не такой случай - - Шаги, как правило, являются смертью для векторных архитектур. Векторные загрузки и магазины не работают таким образом, и много, много работы по загрузке полос в отдельности. В таких случаях вы, как правило, лучше скаляра, хотя некоторые соображения могут позволить вам обыграть скаляр в некоторых случаях с небольшими полями.

То, как вы быстро двигаетесь здесь с векторными свойствами, заключается в том, чтобы сразу найти минимум/максимум для нескольких позиций. Например, если у нас было изображение с плавающей запятой RGBA, тогда найдите min/max для r, g, b, a все одновременно и верните четыре минуты и четыре максимума в конце. Это не намного быстрее, чем код, который у вас есть, но вы получаете больше работы из этого - предполагая, что работа необходима.

Другим методом является сохранение сокращенной копии данных и запуск фильтра по вариантам уменьшенного размера по мере необходимости. Это будет использовать больше памяти, но с коэффициентом два прореживания, он меньше, чем в два раза (1/3 для 2D, c.f. mipmaps).И здесь это полезно, если вы намерены сделать это много.