我需要提高以下代码的性能(Intel Ivy Bridge,x64):
unsigned int delta;
unsigned int a[100];
unsigned int b[100];
...
double sum = 0;
for(int i = 0; i < 100; i++)
sum += (double)b[i]/a[i];
bool possible = delta >= sum;
瓶颈确实是double,并使执行时间增加3倍.a[index]将是从0到500米的任何东西.b[index]将从0到500.
问:在对这段代码的两次调用之间如何修改数组a和b？
在每次调用时,唯一的区别是a[index]++;,其中0 <= index < 100b为总是相同的.三角洲也不会改变.
由于结果与另一个数字进行比较并存储为布尔值,我绝对需要尽可能高的精度.这就是为什么我使用双倍而不是浮动的原因.如你所知,即使1/1b的差异也会返回错误的值,因为结果是布尔值!
小智 8
一件事:
将英特尔组装硬编码到您的程序中将使其不那么便携,更脆弱,更不安全,并且通常是恐怖的工作.这是一项需要避免的任务,除非您需要从裸机中获得最后一盎司的性能,例如编写内核级代码(驱动程序和调度程序).这可能不适合它.
第二件事:
除非你喜欢上帝,否则你可能无法编写比现有代码更快的程序集.C++包含深刻的魔力,许多例程操作都编译成违反直觉的优化,这些优化比天真的解决方案更有效.
三件事:
装配不是你的问题. double表示双精度浮点数,浮点运算通常比整数运算更昂贵,并且这个瓶颈是计算所固有的.
通过考虑算法更改可以最好地解决问题.
您声明更新数组的简单解决方案(减去旧数据并添加新值)是不可接受的,因为精度至关重要.该解决方案具有时间复杂性O(u),其中u是阵列更新的数量和空间复杂性O(1).
到目前为止,所有解决方案都依赖于对整个数组进行重新求和,而每次迭代中只有一个条目发生变化.这是时间复杂性O(un)和空间复杂性O(1).
但显而易见的解决方案是仅重新排列变化的阵列的"部分"!当一个元素在你的数组中更新时,只有一半的数组发生了变化,只有一半的变化了,只有那一半的一半变化了...
我的解决方案是保留每个子阵列总和的完整树.在每次更新时,我都会从叶子向上传播更改的总和,重新使用之前在未更改的子树上完成的子阵列的所有求和.这是以空间复杂O(u log n)性为代价的时间复杂O(n)性.
#include <stdio.h>
#include <time.h>
/**
* Controlling variables.
*/
#ifndef REPEAT
#define REPEAT 2260
#endif
#ifndef NAIVE
#define NAIVE 0
#endif
#ifndef PROBLEM
#define PROBLEM (1<<7)
#endif
#ifndef PRINT_PROGRESS
#define PRINT_PROGRESS 1
#endif
/**
* Initialize the workspace, returning the initial sum.
*/
static double speedyInit(unsigned* a, unsigned* b, unsigned n, double* w, unsigned N){
unsigned i, j;
double* adbl = w+ N;
double* bdbl = w+2*N;
/**
* We initialize the workspace with the correct values out to index i
* and zero(-producing) values from index n to N.
*/
for(i=0;i<n;i++){
adbl[i] = a[i];
bdbl[i] = b[i];
w[i] = bdbl[i]/adbl[i];
}
for(;i<N;i++){
adbl[i] = 1.0;
bdbl[i] = 0.0;
w[i] = 0.0;
}
/**
* We in-place and bottom-up construct the "tree" of sums.
*/
for(i=j=N;i>1;i-=2){/* First-level sums */
w[--j] = w[i-2] + w[i-1];
}
while(--j){/* Subsequent sum levels */
w[j] = w[2*j] + w[2*j+1];
}
/**
* We return the overall sum, found in w[1].
*/
return w[1];
}
/**
* Performs the "A" array update efficiently, returning the new sum.
*/
static double speedyUpd(unsigned* a, double* w, unsigned N, unsigned i){
unsigned p;
double v0, v1;
double* adbl = w+ N;
double* bdbl = w+2*N;
/**
* We increment the two "a" arrays.
*
* NOTE: A double's precision is great enough to losslessly store
* 32-bit unsigned values.
*/
a [i]++;
adbl[i]++;
/**
* We compute the new value at index i and, somewhat wastefully, its "buddy"
* value at index i^1.
*/
v0 = bdbl[i ]/adbl[i ];
v1 = bdbl[i^1]/adbl[i^1];
/**
* We iteratively propagate the v0+v1 sum "up" the top of the "tree" in log-time.
*
* On each iteration we insert the sum v0+v1 at index p, then set v0 to the
* value at index p and v1 to the value of its "buddy", index p^1. The parent
* index of p is then computed and stored in p.
*/
p = (N>>1) + (i>>1);
while(p){
v0 = w[p ] = v0+v1;
v1 = w[p^1];
p >>= 1;
}
/**
* We return the overall sum, found in w[1].
*/
return w[1];
}
/**
* Performs the "A" array update inefficiently, returning the new sum.
*/
static double slowyUpd(unsigned* a, double* w, unsigned N, unsigned i){
double sum = 0;
double* adbl = w+ N;
double* bdbl = w+2*N;
a [i]++;
adbl[i]++;
for(i=0; i<N; i++){
sum += bdbl[i]/adbl[i];
}
return sum;
}
/**
* Requires N a power of two bigger than one.
* Requires n <= N.
* Requires workspace w of 3*N doubles.
*/
double speedy(unsigned* a, unsigned* b, unsigned n, double* w, unsigned N){
int i = n, cond = 1;
double sum;
double delta = 0;
sum = speedyInit(a, b, n, w, N);
while(cond){
/* Do whatever */
/* ... */
/* Set i. */
i = i-1;
/* ... */
#if NAIVE
sum = slowyUpd(a, w, N, rand()%n);
#else
sum = speedyUpd(a, w, N, rand()%n);
#endif
/* ... */
int possible = delta >= sum;
/* ... */
cond = i > 0;
}
return sum;
}
/**
* Main. Gives example.
*/
int main(void){
const unsigned n=PROBLEM, N=PROBLEM;
unsigned a[n], b[n];
double w[3*N];
unsigned i, j;
double dummy = 0;
for(i=0;i<n;i++){
a[i] = 1;
b[i] = i;
}
speedy(a, b, n, w, N);/* Dummy */
clock_t clk = -clock();
for(i=0;i<REPEAT;i++){
dummy += speedy(a, b, n, w, N);
#if PRINT_PROGRESS
putchar('.');
fflush(stdout);
#endif
}
clk += clock();
#if PRINT_PROGRESS
putchar('\n');
#endif
printf("dummy = %f, average time %.9f\n", dummy, clk/((double)CLOCKS_PER_SEC*REPEAT));
}
假设你把它放在名为upd_avg.c的文件中,命令
gcc -O3 upd_avg.c -o upd_avg -DPRINT_PROGRESS=0 -DNAIVE=0 -DREPEAT=2260 -DPROBLEM=128
gcc -O3 upd_avg.c -o upd_avg -DPRINT_PROGRESS=0 -DNAIVE=1 -DREPEAT=2260 -DPROBLEM=128
将分别编译我的O(u log n)算法和其他人的天真O(un)算法.
对于u相同的情况n,差异与日(或mergesort与bubblesort)一样清晰:
| Average time/run (s)
Size | -DNAIVE=0 | -DNAIVE=1
_________________|_____________________________________
-DPROBLEM=2 | 0.000000094 | 0.000000071
-DPROBLEM=4 | 0.000000196 | 0.000000180
-DPROBLEM=8 | 0.000000482 | 0.000000809
-DPROBLEM=16 | 0.000000989 | 0.000002556
... | ... | ...
-DPROBLEM=128 | 0.000007623 | 0.000150181
-DPROBLEM=256 | 0.000016713 | 0.000590156
-DPROBLEM=512 | 0.000037765 | 0.002338671
-DPROBLEM=1024 | 0.000077752 | 0.009324281
-DPROBLEM=2048 | 0.000167924 | 0.037225660
-DPROBLEM=4096 | 0.000343608 | 0.146875721 (*)
... | ... | ...
-DPROBLEM=65536 | 0.007426288 | 21.264978500 (**)
... | ... | We haaaveee liiiffffttttooofffff!!!!!!!
(*) -DREPEAT=226而不是2260.(**)-DREPEAT=2 而不是2260,CPU风扇速度加倍.
我的speedy()函数接受unsigned int数组a和b任何尺寸的n >= 2.然而,它也需要的工作区的存储器的分配,大小的3*N,其中N必须是2的幂,优选等于n四舍五入到2的下一个更高的功率.
该函数speedyInit()设置sums树,从而计算初始值,它位于工作空间的根部,定义w[1]为实现简单的元素.
该函数speedyUpd()是实现我的对数时间和传播的函数.while它内部的循环优雅地实现从树叶向上走树.它由启用-DNAIVE=0.
该功能slowyUpd()是天真的实现.它是由启用的-DNAIVE=1,因为它很慢而被命名.
注意:在尝试对我的代码进行基准测试时,我发现GCC的折叠和DCE在删除无效的代码或仅运行一次函数方面表现出色.
NB我发现精确度至关重要有些奇怪,但是你会依次添加可能大不相同的数量,而不会看到Kahan的求和算法.
| 归档时间： |
|
| 查看次数： |
272 次 |
| 最近记录： |