我想问一下,是否有比逐个迭代所有值并将其除以32768更快的方法来进行音频转换。
void CAudioDataItem::Convert(const vector<int>&uIntegers, vector<double> &uDoubles)
{
for ( int i = 0; i <=uIntegers.size()-1;i++)
{
uDoubles[i] = uIntegers[i] / 32768.0;
}
}
我的方法很好,但可能更快。然而,我没有找到任何加快速度的方法。
谢谢你的帮助!
如果您的数组足够大,那么可能值得将其并行化以用于循环。OpenMP的并行for语句就是我要使用的
然后函数为:
void CAudioDataItem::Convert(const vector<int>&uIntegers, vector<double> &uDoubles)
{
#pragma omp parallel for
for (int i = 0; i < uIntegers.size(); i++)
{
uDoubles[i] = uIntegers[i] / 32768.0;
}
}
对于gcc,在编译要使用的pragma
时需要传递-fopenmp
,在MSVC上它是/openmp
。由于生成线程有明显的开销,所以如果处理大型数组YMMV,这只会更快。
为了获得最大速度,您希望每次循环迭代转换一个以上的值。最简单的方法是使用SIMD。以下是使用SSE2:的大致方法
void CAudioDataItem::Convert(const vector<int>&uIntegers, vector<double> &uDoubles)
{
__m128d scale = _mm_set_pd( 1.0 / 32768.0, 1.0 / 32768.0 );
int i = 0;
for ( ; i < uIntegers.size() - 3; i += 4)
{
__m128i x = _mm_loadu_si128(&uIntegers[i]);
__m128i y = _mm_shuffle_epi32(x, _MM_SHUFFLE(2,3,0,0) );
__m128d dx = _mm_cvtepi32_pd(x);
__m128d dy = _mm_cvtepi32_pd(y);
dx = _mm_mul_pd(dx, scale);
dy = _mm_mul_pd(dy, scale);
_mm_storeu_pd(dx, &uDoubles[i]);
_mm_storeu_pd(dy, &uDoubles[i + 2]);
}
// Finish off the last 0-3 elements the slow way
for ( ; i < uIntegers.size(); i ++)
{
uDoubles[i] = uIntegers[i] / 32768.0;
}
}
我们每次循环迭代处理四个整数。由于寄存器中只能容纳两个double,因此存在一些重复的工作,但除非阵列很小,否则额外的展开将有助于提高性能。
将数据类型更改为较小的数据类型(比如short和float)也应该有助于提高性能,因为它们减少了内存带宽,并且可以在SSE寄存器中容纳四个float。对于音频数据,您不需要双精度。
请注意,我使用了未对齐的加载和存储。如果数据真的对齐了(默认情况下不会对齐,而且很难在std::vector中对齐),那么对齐的速度会稍微快一点。
您的函数具有高度并行性。在现代英特尔CPU上,有三种独立的并行方式:指令级并行(ILP)、线程级并行(TLP)和SIMD我能够用这三个来大大提升你的功能但结果取决于编译器。使用GCC的提升要少得多,因为它已经对函数进行了矢量化。请参阅下面的数字表。
然而,函数中的主要限制因素是它的时间复杂度仅为O(n),因此当您运行的阵列大小跨越每个缓存级别边界时,效率会急剧下降。例如,如果你看一下大型密集矩阵乘法(GEMM),它是一个O(n^3)运算,所以如果一个人做得对(例如使用循环平铺),缓存层次结构就不是问题:即使对于非常大的矩阵,你也可以接近最大flop/s(这似乎表明GEMM是英特尔设计CPU时考虑的问题之一)。在您的情况下,解决此问题的方法是在执行更复杂的操作(例如,O(n^2))之后/之前,找到一种方法在一级缓存块上执行功能,然后移动到另一个一级块。我当然不知道你在做什么,所以我不能那样做。
ILP部分由CPU硬件为您完成。然而,经常携带的循环依赖关系限制了ILP,因此进行循环展开以充分利用ILP通常会有所帮助。对于TLP,我使用OpenMP,对于SIMD,我使用AVX(但是下面的代码也适用于SSE)。我使用了32字节对齐的内存,并确保数组是8的倍数,这样就不需要清理了。
以下是Visual Studio 2012 64位AVX和OpenMP(显然是发布模式)的结果SandyBridge EP 4核(8个HW线程)@3.6 GHz。变量n
是项目数。我也重复了好几次函数,所以总时间包括了这一点。函数convert_vec4_unroll2_openmp
给出了除L1区域之外的最佳结果。您还可以清楚地看到,每次移动到新的缓存级别时,效率都会显著下降,但即使是主内存,效率也会更好。
l1 chache, n 2752, repeat 300000
covert time 1.34, error 0.000000
convert_vec4 time 0.16, error 0.000000
convert_vec4_unroll2 time 0.16, error 0.000000
convert_vec4_unroll2_openmp time 0.31, error 0.000000
l2 chache, n 21856, repeat 30000
covert time 1.14, error 0.000000
convert_vec4 time 0.24, error 0.000000
convert_vec4_unroll2 time 0.24, error 0.000000
convert_vec4_unroll2_openmp time 0.12, error 0.000000
l3 chache, n 699072, repeat 1000
covert time 1.23, error 0.000000
convert_vec4 time 0.44, error 0.000000
convert_vec4_unroll2 time 0.45, error 0.000000
convert_vec4_unroll2_openmp time 0.14, error 0.000000
main memory , n 8738144, repeat 100
covert time 1.56, error 0.000000
convert_vec4 time 0.95, error 0.000000
convert_vec4_unroll2 time 0.89, error 0.000000
convert_vec4_unroll2_openmp time 0.51, error 0.000000
g++ foo.cpp -mavx -fopenmp -ffast-math -O3
在i5-3317(常青藤桥)上的结果@2.4 GHz 2核(4个HW线程)。GCC似乎将其矢量化,唯一的好处来自OpenMP(然而,它在L1区域给出了更糟糕的结果)。
l1 chache, n 2752, repeat 300000
covert time 0.26, error 0.000000
convert_vec4 time 0.22, error 0.000000
convert_vec4_unroll2 time 0.21, error 0.000000
convert_vec4_unroll2_openmp time 0.46, error 0.000000
l2 chache, n 21856, repeat 30000
covert time 0.28, error 0.000000
convert_vec4 time 0.27, error 0.000000
convert_vec4_unroll2 time 0.27, error 0.000000
convert_vec4_unroll2_openmp time 0.20, error 0.000000
l3 chache, n 699072, repeat 1000
covert time 0.80, error 0.000000
convert_vec4 time 0.80, error 0.000000
convert_vec4_unroll2 time 0.80, error 0.000000
convert_vec4_unroll2_openmp time 0.83, error 0.000000
main memory chache, n 8738144, repeat 100
covert time 1.10, error 0.000000
convert_vec4 time 1.10, error 0.000000
convert_vec4_unroll2 time 1.10, error 0.000000
convert_vec4_unroll2_openmp time 1.00, error 0.000000
这是代码。我使用vectorclasshttp://www.agner.org/optimize/vectorclass.zip进行SIMD。这将使用AVX一次写入4个双精度,或SSE一次写入2个双精度。
#include <stdlib.h>
#include <stdio.h>
#include <omp.h>
#include "vectorclass.h"
void convert(const int *uIntegers, double *uDoubles, const int n) {
for ( int i = 0; i<n; i++) {
uDoubles[i] = uIntegers[i] / 32768.0;
}
}
void convert_vec4(const int *uIntegers, double *uDoubles, const int n) {
Vec4d div = 1.0/32768;
for ( int i = 0; i<n; i+=4) {
Vec4i u4i = Vec4i().load(&uIntegers[i]);
Vec4d u4d = to_double(u4i);
u4d*=div;
u4d.store(&uDoubles[i]);
}
}
void convert_vec4_unroll2(const int *uIntegers, double *uDoubles, const int n) {
Vec4d div = 1.0/32768;
for ( int i = 0; i<n; i+=8) {
Vec4i u4i_v1 = Vec4i().load(&uIntegers[i]);
Vec4d u4d_v1 = to_double(u4i_v1);
u4d_v1*=div;
u4d_v1.store(&uDoubles[i]);
Vec4i u4i_v2 = Vec4i().load(&uIntegers[i+4]);
Vec4d u4d_v2 = to_double(u4i_v2);
u4d_v2*=div;
u4d_v2.store(&uDoubles[i+4]);
}
}
void convert_vec4_openmp(const int *uIntegers, double *uDoubles, const int n) {
#pragma omp parallel for
for ( int i = 0; i<n; i+=4) {
Vec4i u4i = Vec4i().load(&uIntegers[i]);
Vec4d u4d = to_double(u4i);
u4d/=32768.0;
u4d.store(&uDoubles[i]);
}
}
void convert_vec4_unroll2_openmp(const int *uIntegers, double *uDoubles, const int n) {
Vec4d div = 1.0/32768;
#pragma omp parallel for
for ( int i = 0; i<n; i+=8) {
Vec4i u4i_v1 = Vec4i().load(&uIntegers[i]);
Vec4d u4d_v1 = to_double(u4i_v1);
u4d_v1*=div;
u4d_v1.store(&uDoubles[i]);
Vec4i u4i_v2 = Vec4i().load(&uIntegers[i+4]);
Vec4d u4d_v2 = to_double(u4i_v2);
u4d_v2*=div;
u4d_v2.store(&uDoubles[i+4]);
}
}
double compare(double *a, double *b, const int n) {
double diff = 0.0;
for(int i=0; i<n; i++) {
double tmp = a[i] - b[i];
//printf("%d %f %f n", i, a[i], b[i]);
if(tmp<0) tmp*=-1;
diff += tmp;
}
return diff;
}
void loop(const int n, const int repeat, const int ifunc) {
void (*fp[4])(const int *uIntegers, double *uDoubles, const int n);
int *a = (int*)_mm_malloc(sizeof(int)* n, 32);
double *b1_cmp = (double*)_mm_malloc(sizeof(double)*n, 32);
double *b1 = (double*)_mm_malloc(sizeof(double)*n, 32);
double dtime;
const char *fp_str[] = {
"covert",
"convert_vec4",
"convert_vec4_unroll2",
"convert_vec4_unroll2_openmp",
};
for(int i=0; i<n; i++) {
a[i] = rand()*RAND_MAX;
}
fp[0] = convert;
fp[1] = convert_vec4;
fp[2] = convert_vec4_unroll2;
fp[3] = convert_vec4_unroll2_openmp;
convert(a, b1_cmp, n);
dtime = omp_get_wtime();
for(int i=0; i<repeat; i++) {
fp[ifunc](a, b1, n);
}
dtime = omp_get_wtime() - dtime;
printf("t%s time %.2f, error %fn", fp_str[ifunc], dtime, compare(b1_cmp,b1,n));
_mm_free(a);
_mm_free(b1_cmp);
_mm_free(b1);
}
int main() {
double dtime;
int l1 = (32*1024)/(sizeof(int) + sizeof(double));
int l2 = (256*1024)/(sizeof(int) + sizeof(double));
int l3 = (8*1024*1024)/(sizeof(int) + sizeof(double));
int lx = (100*1024*1024)/(sizeof(int) + sizeof(double));
int n[] = {l1, l2, l3, lx};
int repeat[] = {300000, 30000, 1000, 100};
const char *cache_str[] = {"l1", "l2", "l3", "main memory"};
for(int c=0; c<4; c++ ) {
int lda = ((n[c]+7) & -8); //make sure array is a multiple of 8
printf("%s chache, n %dn", cache_str[c], lda);
for(int i=0; i<4; i++) {
loop(lda, repeat[c], i);
} printf("n");
}
}
最后,任何读过这篇文章并想提醒我我的代码看起来更像C而不是C++的人,请在决定发表评论之前先阅读这篇文章http://www.stroustrup.com/sibling_rivalry.pdf
您也可以尝试:
uDoubles[i] = ldexp((double)uIntegers[i], -15);
编辑:有关使用SSE内部函数的版本,请参阅上面Adam的回答。比我在这里更好。。。
为了使其更加有用,让我们在这里查看编译器生成的代码。我使用的是gcc 4.8.0,是的,值得检查您的特定编译器(版本),因为gcc 4.4、4.8、clang 3.2或英特尔的icc的输出存在相当大的差异。
您的原件,使用g++ -O8 -msse4.2 ...
翻译成以下循环:
.L2:
cvtsi2sd (%rcx,%rax,4), %xmm0
mulsd %xmm1, %xmm0
addl $1, %edx
movsd %xmm0, (%rsi,%rax,8)
movslq %edx, %rax
cmpq %rdi, %rax
jbe .L2
其中CCD_ 8保持CCD_。
另一方面,使用g++ -msse4.2 -O8 -funroll-loops ...
,为循环创建的代码发生了显著变化:
[ ... ]
leaq -1(%rax), %rdi
movq %rdi, %r8
andl $7, %r8d
je .L3
[ ... insert a duff's device here, up to 6 * 2 conversions ... ]
jmp .L3
.p2align 4,,10
.p2align 3
.L39:
leaq 2(%rsi), %r11
cvtsi2sd (%rdx,%r10,4), %xmm9
mulsd %xmm0, %xmm9
leaq 5(%rsi), %r9
leaq 3(%rsi), %rax
leaq 4(%rsi), %r8
cvtsi2sd (%rdx,%r11,4), %xmm10
mulsd %xmm0, %xmm10
cvtsi2sd (%rdx,%rax,4), %xmm11
cvtsi2sd (%rdx,%r8,4), %xmm12
cvtsi2sd (%rdx,%r9,4), %xmm13
movsd %xmm9, (%rcx,%r10,8)
leaq 6(%rsi), %r10
mulsd %xmm0, %xmm11
mulsd %xmm0, %xmm12
movsd %xmm10, (%rcx,%r11,8)
leaq 7(%rsi), %r11
mulsd %xmm0, %xmm13
cvtsi2sd (%rdx,%r10,4), %xmm14
mulsd %xmm0, %xmm14
cvtsi2sd (%rdx,%r11,4), %xmm15
mulsd %xmm0, %xmm15
movsd %xmm11, (%rcx,%rax,8)
movsd %xmm12, (%rcx,%r8,8)
movsd %xmm13, (%rcx,%r9,8)
leaq 8(%rsi), %r9
movsd %xmm14, (%rcx,%r10,8)
movsd %xmm15, (%rcx,%r11,8)
movq %r9, %rsi
.L3:
cvtsi2sd (%rdx,%r9,4), %xmm8
mulsd %xmm0, %xmm8
leaq 1(%rsi), %r10
cmpq %rdi, %r10
movsd %xmm8, (%rcx,%r9,8)
jbe .L39
[ ... out ... ]
因此,它阻止了操作,但仍然一次转换一个值。
如果您更改原始循环,以便在每次迭代中对几个元素进行操作:
size_t i;
for (i = 0; i < uIntegers.size() - 3; i += 4)
{
uDoubles[i] = uIntegers[i] / 32768.0;
uDoubles[i+1] = uIntegers[i+1] / 32768.0;
uDoubles[i+2] = uIntegers[i+2] / 32768.0;
uDoubles[i+3] = uIntegers[i+3] / 32768.0;
}
for (; i < uIntegers.size(); i++)
uDoubles[i] = uIntegers[i] / 32768.0;
编译器gcc -msse4.2 -O8 ...
(即即使没有请求展开)识别出使用CVTDQ2PD
/MULPD
的可能性,循环的核心变成:
.p2align 4,,10
.p2align 3
.L4:
movdqu (%rcx), %xmm0
addq $16, %rcx
cvtdq2pd %xmm0, %xmm1
pshufd $238, %xmm0, %xmm0
mulpd %xmm2, %xmm1
cvtdq2pd %xmm0, %xmm0
mulpd %xmm2, %xmm0
movlpd %xmm1, (%rdx,%rax,8)
movhpd %xmm1, 8(%rdx,%rax,8)
movlpd %xmm0, 16(%rdx,%rax,8)
movhpd %xmm0, 24(%rdx,%rax,8)
addq $4, %rax
cmpq %r8, %rax
jb .L4
cmpq %rdi, %rax
jae .L29
[ ... duff's device style for the "tail" ... ]
.L29:
rep ret
也就是说,现在编译器认识到有机会在每个SSE寄存器中放入两个double
,并进行并行乘法/转换。这与Adam的SSE内部版本将生成的代码非常接近。
总的来说,代码(我只展示了大约1/6)比"直接"内部函数要复杂得多,因为如上所述,编译器试图在循环中预写/附加未对齐/未阻塞多个"头"one_answers"尾"。这在很大程度上取决于矢量的平均/预期大小,这是否有益;对于"通用"情况(向量是"最内部"循环处理的块大小的两倍以上),它会有所帮助。
这次演习的结果在很大程度上是。。。如果你强迫(通过编译器选项/优化)或暗示(通过稍微重新排列代码)你的编译器做正确的事情,那么对于这种特定类型的复制/转换循环,它会产生不会落后于手工编写的内部函数的代码。
最后的实验。。。生成代码:
static double c(int x) { return x / 32768.0; }
void Convert(const std::vector<int>& uIntegers, std::vector<double>& uDoubles)
{
std::transform(uIntegers.begin(), uIntegers.end(), uDoubles.begin(), c);
}
并且(为了获得最好的汇编输出,这次使用带有gcc -O8 -msse4.2 ...
的gcc 4.4)生成的汇编核心循环(同样,有一个前置/后置位)变为:
.p2align 4,,10
.p2align 3
.L8:
movdqu (%r9,%rax), %xmm0
addq $1, %rcx
cvtdq2pd %xmm0, %xmm1
pshufd $238, %xmm0, %xmm0
mulpd %xmm2, %xmm1
cvtdq2pd %xmm0, %xmm0
mulpd %xmm2, %xmm0
movapd %xmm1, (%rsi,%rax,2)
movapd %xmm0, 16(%rsi,%rax,2)
addq $16, %rax
cmpq %rcx, %rdi
ja .L8
cmpq %rbx, %rbp
leaq (%r11,%rbx,4), %r11
leaq (%rdx,%rbx,8), %rdx
je .L10
[ ... ]
.L10:
[ ... ]
ret
有了这些,我们能学到什么?如果你想使用C++,真的要使用C++;-)
让我尝试另一种方法:
如果从汇编指令的角度来看,乘法是非常好的,那么这应该保证它会被相乘。
void CAudioDataItem::Convert(const vector<int>&uIntegers, vector<double> &uDoubles)
{
for ( int i = 0; i <=uIntegers.size()-1;i++)
{
uDoubles[i] = uIntegers[i] * 0.000030517578125;
}
}