我想问一下是否有比逐个遍历所有值并将它们除以 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;
}
}
我的方法效果很好，但可能会更快。但是我没有找到任何方法来加快它。
感谢您的帮助！
如果您的数组足够大，那么并行化这个 for 循环可能是值得的。我会使用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 时，您需要-fopenmp在编译pragma以供使用时传递，在 MSVC 上它是/openmp. 由于生成线程具有明显的开销，因此只有在处理大型数组 YMMV 时才会更快。
您的函数是高度可并行化的。在现代 Intel CPU 上，有三种独立的并行方式：指令级并行 (ILP)、线程级并行 (TLP) 和 SIMD。 我能够使用所有这三个来大大提升您的功能。 结果取决于编译器。使用 GCC 的提升要少得多，因为它已经对函数进行了矢量化。请参阅下面的数字表。
但是，函数的主要限制因素是它的时间复杂度仅为 O(n)，因此当您运行的数组的大小跨越每个缓存级别边界时，效率会急剧下降. 例如，如果您看一下大型密集矩阵乘法 (GEMM)，它是一个 O(n^3) 操作，因此如果一个操作正确（使用例如循环平铺），则缓存层次结构不是问题：您可以接近最大值即使对于非常大的矩阵也有 flops/s（这似乎表明 GEMM 是英特尔在设计 CPU 时考虑的想法之一）。在您的情况下解决此问题的方法是在您执行更复杂的操作之后/之前找到一种方法在 L1 缓存块上执行您的功能（例如，作为 O(n^2)）然后移动到另一个L1 块。当然，我不知道你在做什么，所以我不能那样做。
ILP 部分由 CPU 硬件为您完成。然而，经常携带的循环依赖限制了 ILP，因此它通常有助于进行循环展开以充分利用 ILP。对于 TLP，我使用 OpenMP，对于 SIMD，我使用 AVX（但是下面的代码也适用于 SSE）。我使用了 32 字节对齐的内存，并确保数组是 8 的倍数，因此不需要清理。
这是 Visual Studio 2012 64bit 的结果，带有 AVX 和 OpenMP（显然是发布模式）SandyBridge EP 4 核（8 个硬件线程）@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
在i5-3317 g++ foo.cpp -mavx -fopenmp -ffast-math -O3（常春藤桥）@ 2.4 GHz 2 核（4 个硬件线程）上的结果。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
这是代码。我使用vectorclass http://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 %f\n", 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 %d\n", 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
为了获得最大速度，您希望每次循环迭代转换多个值。最简单的方法是使用 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;
}
}
我们每次循环迭代处理四个整数。由于我们只能在寄存器中放置两个双精度数，因此存在一些重复的工作，但是除非数组很小，否则额外的展开将有助于提高性能。
将数据类型更改为更小的数据类型（比如 short 和 float）也应该有助于提高性能，因为它们会减少内存带宽，并且您可以在 SSE 寄存器中放置四个浮点数。对于音频数据，您不需要双精度。
请注意，我使用了未对齐的加载和存储。如果数据实际上是对齐的，对齐的会稍微快一些（默认情况下不会对齐，并且很难在 std::vector 内对齐）。
你也可以试试：
uDoubles[i] = ldexp((double)uIntegers[i], -15);
编辑：有关使用 SSE 内在函数的版本，请参阅上面的 Adam 的答案。比我这里的好...
为了使它更有用，让我们在这里看看编译器生成的代码。我正在使用 gcc 4.8.0，是的，值得检查您的特定编译器（版本），因为 gcc 4.4、4.8、clang 3.2 或 Intel 的 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
where%xmm1成立1.0/32768.0，因此编译器会自动将除法转换为反向乘法。
另一方面，使用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
即现在编译器认识到有机会在double每个 SSE 寄存器中放置两个，并进行并行乘法/转换。这与 Adam 的 SSE 内在函数版本将生成的代码非常接近。
总代码（我只展示了其中的 1/6）比“直接”内在函数复杂得多，因为如前所述，编译器尝试预先添加/附加未对齐/非块-循环的多个“头”和“尾”。这在很大程度上取决于向量的平均/预期大小，这是否有益；对于“通用”情况（向量是“最内层”循环处理的块大小的两倍以上），它会有所帮助。
这个练习的结果主要是......如果你强制（通过编译器选项/优化）或提示（通过稍微重新排列代码）你的编译器做正确的事情，那么对于这种特定类型的复制/转换循环，它提供的代码不会落后于手写的内在函数。
最后的实验......制作代码：
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 4.4 with gcc -O8 -msse4.2 ...）生成的程序集核心循环（同样，有一个前/后位）变为：
.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;
}
}