You can use SSE2/SSE3 minps / pminsd or relevant instruction set for your processor/architectoure since it is supported directly in GCC / MASM / TASM (In case MASM or TASM is not supported such SSE2/SSE3 instruction set there are also the .inc files to create macros simulating instruction sets on the web for MASM), create .OBJ file by your favorite linker then link it as usual and use in you favorite IDE. You will get from 4x to 16x performance boost compared to the traditional "classic" algorithm. It depend on data size (old compilers treats double not in IEEE format, bout like float in several configurations, on 16x systems, particularly, double means 32 bit data structure, not 64 bit data structure, in modern languages it is correlated to "double" and "long double" data structures, respectively)
The idea is simple: if you have k elements, [k=4n+p, 4>p=>0], complete it with n-p elements or just load last 4 doubles resetting to 0 last p elements, so you can fast evaluate n candidates. evaluate candidates n times comparing to the accumulator, you will get a minimum.
If your processor supports SSE5 or is a brand new, most likely you also will be using one of the HD instructions, which really handy, because it can find maximum (not minimum yet) in array of double values.
Sample of using SSE to calculate peak values of a float array:
#include <xmmintrin.h>
double min(double* array, int size) {
// returns the minimum value of array
int i;
double val = array[0];
for (i = 1; i < size; ++i) {
if (val > array[i]) {
val = array[i];
}
}
return val;
}
#define ARRAY_SIZE 16
float m_fArray[ARRAY_SIZE];
void x86_sse_find_peaks(float *buf, unsigned nframes, float *min, float *max)
{
__m128 current_max, current_min, work;
// Load max and min values into all four slots of the XMM registers
current_min = _mm_set1_ps(*min);
current_max = _mm_set1_ps(*max);
// Work input until "buf" reaches 16 byte alignment
while ( ((unsigned long)buf) % 16 != 0 && nframes > 0) {
// Load the next float into the work buffer
work = _mm_set1_ps(*buf);
current_min = _mm_min_ps(current_min, work);
current_max = _mm_max_ps(current_max, work);
buf++;
nframes--;
}
// use 64 byte prefetch for quadruple quads
while (nframes >= 16) {
//__builtin_prefetch(buf+64,0,0); // for GCC 4.3.2+
work = _mm_load_ps(buf);
current_min = _mm_min_ps(current_min, work);
current_max = _mm_max_ps(current_max, work);
buf+=4;
work = _mm_load_ps(buf);
current_min = _mm_min_ps(current_min, work);
current_max = _mm_max_ps(current_max, work);
buf+=4;
work = _mm_load_ps(buf);
current_min = _mm_min_ps(current_min, work);
current_max = _mm_max_ps(current_max, work);
buf+=4;
work = _mm_load_ps(buf);
current_min = _mm_min_ps(current_min, work);
current_max = _mm_max_ps(current_max, work);
buf+=4;
nframes-=16;
}
// work through aligned buffers
while (nframes >= 4) {
work = _mm_load_ps(buf);
current_min = _mm_min_ps(current_min, work);
current_max = _mm_max_ps(current_max, work);
buf+=4;
nframes-=4;
}
// work through the rest < 4 samples
while ( nframes > 0) {
// Load the next float into the work buffer
work = _mm_set1_ps(*buf);
current_min = _mm_min_ps(current_min, work);
current_max = _mm_max_ps(current_max, work);
buf++;
nframes--;
}
// Find min & max value in current_max through shuffle tricks
work = current_min;
work = _mm_shuffle_ps(work, work, _MM_SHUFFLE(2, 3, 0, 1));
work = _mm_min_ps (work, current_min);
current_min = work;
work = _mm_shuffle_ps(work, work, _MM_SHUFFLE(1, 0, 3, 2));
work = _mm_min_ps (work, current_min);
_mm_store_ss(min, work);
work = current_max;
work = _mm_shuffle_ps(work, work, _MM_SHUFFLE(2, 3, 0, 1));
work = _mm_max_ps (work, current_max);
current_max = work;
work = _mm_shuffle_ps(work, work, _MM_SHUFFLE(1, 0, 3, 2));
work = _mm_max_ps (work, current_max);
_mm_store_ss(max, work);
}
int _tmain(int argc, _TCHAR* argv[])
{
float min = FLT_MAX;
float max = FLT_MIN;
m_fArray[0] = 0;
m_fArray[1] = 1;
m_fArray[2] = 2;
m_fArray[3] = 3;
m_fArray[4] = 4;
m_fArray[5] = 3;
m_fArray[6] = 2;
m_fArray[7] = 1;
m_fArray[8] = -1;
m_fArray[9] = -2;
m_fArray[10] = -3;
m_fArray[11] = -4;
m_fArray[12] = -5;
m_fArray[13] = -6;
m_fArray[14] = -7;
m_fArray[15] = -8;
x86_sse_find_peaks(m_fArray, ARRAY_SIZE, &min, &max);
printf("value = %.2f, max = %.2f\n", min, max); // output is: value = -8.00, max = 4.00
return 0;
}
This needs to be robust
\$\endgroup\$