Здравствуйте, я пытаюсь запустить программу, которая находит ближайшую пару, используя грубую силу с методами кэширования, такими как PDF здесь: Кэширование производительности Стэнфорд
Мой оригинальный код:
float compare_points_BF(int N,point *P){
int i,j;
float distance=0, min_dist=FLT_MAX;
point *p1, *p2;
unsigned long long calc = 0;
for (i=0;i<(N-1);i++){
for (j=i+1;j<N;j++){
if ((distance = (P[i].x - P[j].x) * (P[i].x - P[j].x) +
(P[i].y - P[j].y) * (P[i].y - P[j].y)) < min_dist){
min_dist = distance;
p1 = &P[i];
p2 = &P[j];
}
}
}
return sqrt(min_dist);
}
Эта программа дает примерно эти времена работы:
N 8192 16384 32768 65536 131072 262144 524288 1048576
seconds 0,070 0,280 1,130 5,540 18,080 72,838 295,660 1220,576
0,080 0,330 1,280 5,190 20,290 80,880 326,460 1318,631
Кеш-версия вышеуказанной программы:
float compare_points_BF(register int N, register int B, point *P){
register int i, j, ib, jb, num_blocks = (N + (B-1)) / B;
register point *p1, *p2;
register float distance=0, min_dist=FLT_MAX, regx, regy;
//break array data in N/B blocks, ib is index for i cached block and jb is index for j strided cached block
//each i block is compared with the j block, (which j block is always after the i block)
for (i = 0; i < num_blocks; i++){
for (j = i; j < num_blocks; j++){
//reads the moving frame block to compare with the i cached block
for (jb = j * B; jb < ( ((j+1)*B) < N ? ((j+1)*B) : N); jb++){
//avoid float comparisons that occur when i block = j block
//Register Allocated
regx = P[jb].x;
regy = P[jb].y;
for (i == j ? (ib = jb + 1) : (ib = i * B); ib < ( ((i+1)*B) < N ? ((i+1)*B) : N); ib++){
//calculate distance of current points
if((distance = (P[ib].x - regx) * (P[ib].x - regx) +
(P[ib].y - regy) * (P[ib].y - regy)) < min_dist){
min_dist = distance;
p1 = &P[ib];
p2 = &P[jb];
}
}
}
}
}
return sqrt(min_dist);
}
и некоторые результаты:
Block_size = 256 N = 8192 Run time: 0.090 sec
Block_size = 512 N = 8192 Run time: 0.090 sec
Block_size = 1024 N = 8192 Run time: 0.090 sec
Block_size = 2048 N = 8192 Run time: 0.100 sec
Block_size = 4096 N = 8192 Run time: 0.090 sec
Block_size = 8192 N = 8192 Run time: 0.090 secBlock_size = 256 N = 16384 Run time: 0.357 sec
Block_size = 512 N = 16384 Run time: 0.353 sec
Block_size = 1024 N = 16384 Run time: 0.360 sec
Block_size = 2048 N = 16384 Run time: 0.360 sec
Block_size = 4096 N = 16384 Run time: 0.370 sec
Block_size = 8192 N = 16384 Run time: 0.350 sec
Block_size = 16384 N = 16384 Run time: 0.350 sec
Block_size = 128 N = 32768 Run time: 1.420 sec
Block_size = 256 N = 32768 Run time: 1.420 sec
Block_size = 512 N = 32768 Run time: 1.390 sec
Block_size = 1024 N = 32768 Run time: 1.410 sec
Block_size = 2048 N = 32768 Run time: 1.430 sec
Block_size = 4096 N = 32768 Run time: 1.430 sec
Block_size = 8192 N = 32768 Run time: 1.400 sec
Block_size = 16384 N = 32768 Run time: 1.380 sec
Block_size = 256 N = 65536 Run time: 5.760 sec
Block_size = 512 N = 65536 Run time: 5.790 sec
Block_size = 1024 N = 65536 Run time: 5.720 sec
Block_size = 2048 N = 65536 Run time: 5.720 sec
Block_size = 4096 N = 65536 Run time: 5.720 sec
Block_size = 8192 N = 65536 Run time: 5.530 sec
Block_size = 16384 N = 65536 Run time: 5.550 sec
Block_size = 256 N = 131072 Run time: 22.750 sec
Block_size = 512 N = 131072 Run time: 23.130 sec
Block_size = 1024 N = 131072 Run time: 22.810 sec
Block_size = 2048 N = 131072 Run time: 22.690 sec
Block_size = 4096 N = 131072 Run time: 22.710 sec
Block_size = 8192 N = 131072 Run time: 21.970 sec
Block_size = 16384 N = 131072 Run time: 22.010 sec
Block_size = 256 N = 262144 Run time: 90.220 sec
Block_size = 512 N = 262144 Run time: 92.140 sec
Block_size = 1024 N = 262144 Run time: 91.181 sec
Block_size = 2048 N = 262144 Run time: 90.681 sec
Block_size = 4096 N = 262144 Run time: 90.760 sec
Block_size = 8192 N = 262144 Run time: 87.660 sec
Block_size = 16384 N = 262144 Run time: 87.760 sec
Block_size = 256 N = 524288 Run time: 361.151 sec
Block_size = 512 N = 524288 Run time: 379.521 sec
Block_size = 1024 N = 524288 Run time: 379.801 sec
Из того, что мы видим, время выполнения медленнее, чем у не кешированного кода.
Это связано с оптимизацией компилятора? Это плохой код или это просто из-за алгоритма, который не очень хорошо работает с тайлингом? Я использую VS 2010, скомпилированный с 32-битным исполняемым файлом. Заранее спасибо!
Это интересный случай. Компилятор плохо справился с подъемом инварианта цикла в двух внутренних циклах. А именно, два внутренних цикла for проверяют следующее условие в каждой итерации:
(j+1)*B) < N ? ((j+1)*B) : N
а также
(i+1)*B) < N ? ((i+1)*B) : N
Расчет и ветвление дорогие; но они на самом деле являются инвариантами цикла для двух внутренних циклов for. После того, как я вручную вытащил их из двух внутренних циклов for, я смог получить оптимизированную кеш-версию, чтобы она работала лучше, чем неоптимизированная (10% при N == 524288, 30% при N = 1048576).
Вот модифицированный код (простое изменение, ищите u1, u2):
//break array data in N/B blocks, ib is index for i cached block and jb is index for j strided cached block
//each i block is compared with the j block, (which j block is always after the i block)
for (i = 0; i < num_blocks; i++){
for (j = i; j < num_blocks; j++){
int u1 = (((j+1)*B) < N ? ((j+1)*B) : N);
int u2 = (((i+1)*B) < N ? ((i+1)*B) : N);
//reads the moving frame block to compare with the i cached block
for (jb = j * B; jb < u1 ; jb++){
//avoid float comparisons that occur when i block = j block
//Register Allocated
regx = P[jb].x;
regy = P[jb].y;
for (i == j ? (ib = jb + 1) : (ib = i * B); ib < u2; ib++){
//calculate distance of current points
if((distance = (P[ib].x - regx) * (P[ib].x - regx) +
(P[ib].y - regy) * (P[ib].y - regy)) < min_dist){
min_dist = distance;
p1 = &P[ib];
p2 = &P[jb];
}
}
}
}
}
Черепица может быть старая концепция, но она все еще очень актуальна сегодня. В исходном фрагменте кода для каждого i вы можете повторно использовать большинство элементов P [j], пока они все еще кэшируются, но только в том случае, если длина внутреннего цикла была достаточно мала, чтобы поместиться там. Фактический размер должен быть определен тем, какой уровень кэша вы хотите использовать для разбиения на листы — L1 обеспечит наилучшую производительность, поскольку он самый быстрый, но так как он также наименьший, вам понадобятся небольшие блоки, а накладные расходы могут быть слишком большими. L2 допускает большие плитки, но незначительно снижает производительность и так далее.
Обратите внимание, что вам не нужно использовать двухмерное разбиение здесь, это не матричное умножение — вы пересекаете один и тот же массив. Вы можете просто разбить внутренний цикл, так как он переполняет кеш, как только вы это сделаете — внешний цикл (i) может выполняться до конца текущего блока кэшированных элементов внутреннего цикла. На самом деле нет никакого смысла в двумерном тайлинге, поскольку никто не собирается повторно использовать элементы внешнего цикла (в отличие от матрицы mul)
Итак, предполагая Point
имеет размер 64 бита, вы можете безопасно разместить 512 таких элементов массива в 32k L1 или 4096 элементов в 256k L2. вам придется пропустить один раз для P [i] в каждом блоке, если i выходит за пределы текущего j-го блока, но это незначительно.
Кстати, это объяснение все еще может быть устаревшим, так как достаточно хороший компилятор может попытаться сделать все это за вас. Хотя это довольно сложно, поэтому я немного скептически отношусь к любому из распространенных, но здесь должно быть легко доказать, что переупорядочение безопасно. Можно, конечно, утверждать, что «достаточно хороший компилятор» — это парадокс, но это не по теме …