如何从具有向量数组位置的数组(内存中)中以最佳方式读取?

Question

我有这样的代码：

const rack::simd::float_4 pos = phase * waveTable.mLength;
const rack::simd::int32_4 pos0 = pos;
const rack::simd::float_4 frac = pos - (rack::simd::float_4)pos0;
rack::simd::float_4 v0;
rack::simd::float_4 v1;
for (int v = 0; v < 4; v++) {
v0[v] = waveTable.mSamples[pos0[v]];
v1[v] = waveTable.mSamples[pos0[v] + 1]; // mSamples size is waveTable.mLength + 1 for interpolation wraparound
}
oversampleBuffer[i] = v0 + (v1 - v0) * frac;

它采用两个样本之间的phase(归一化)和插值(线性插值)，存储在waveTable.mSamples上(每个样本作为单个浮点数)。

位置是rack::simd::float_4内部，基本上是4对齐的浮点数定义为__m128。 经过一些基准测试后，这部分代码需要一些时间(我猜是由于缺少大量缓存)。

是用-march=nocona构建的，所以我可以使用MMX，SSE，SSE2和SSE3。

您将如何优化此代码？谢谢

Answer 1

由于多种原因，您的代码效率不太高。

1. 您正在使用标量代码设置 SIMD 矢量的各个通道。处理器不能完全做到这一点，但编译器假装他们可以。不幸的是，这些编译器实现的解决方法很慢，通常它们通过对内存和返回进行漫游来做到这一点。

2. 通常，应避免写入 2 或 4 个非常小长度的循环。有时编译器展开，你没问题，但有时他们没有，CPU错误地预测了太多的分支。

3. 最后，处理器可以使用单个指令加载 64 位值。您正在从表中加载连续对，可以使用 64 位加载而不是两个 32 位加载。

这是一个固定版本(未经测试)。这假设您正在为 PC 构建，即使用 SSE SIMD。

// Load a vector with rsi[ i0 ], rsi[ i0 + 1 ], rsi[ i1 ], rsi[ i1 + 1 ]
inline __m128 loadFloats( const float* rsi, int i0, int i1 )
{
// Casting load indices to unsigned, otherwise compiler will emit sign extension instructions
__m128d res = _mm_load_sd( (const double*)( rsi + (uint32_t)i0 ) );
res = _mm_loadh_pd( res, (const double*)( rsi + (uint32_t)i1 ) );
return _mm_castpd_ps( res );
}
__m128 interpolate4( const float* waveTableData, uint32_t waveTableLength, __m128 phase )
{
// Convert wave table length into floats.
// Consider doing that outside of the inner loop, and passing the __m128.
const __m128 length = _mm_set1_ps( (float)waveTableLength );
// Compute integer indices, and the fraction
const __m128 pos = _mm_mul_ps( phase, length );
const __m128 posFloor = _mm_floor_ps( pos );    // BTW this one needs SSE 4.1, workarounds are expensive
const __m128 frac = _mm_sub_ps( pos, posFloor );
const __m128i posInt = _mm_cvtps_epi32( posFloor );
// Abuse 64-bit load instructions to load pairs of values from the table.
// If you have AVX2, can use _mm256_i32gather_pd instead, will load all 8 floats with 1 (slow) instruction.
const __m128 s01 = loadFloats( waveTableData, _mm_cvtsi128_si32( posInt ), _mm_extract_epi32( posInt, 1 ) );
const __m128 s23 = loadFloats( waveTableData, _mm_extract_epi32( posInt, 2 ), _mm_extract_epi32( posInt, 3 ) );
// Shuffle into the correct order, i.e. gather even/odd elements from the vectors
const __m128 v0 = _mm_shuffle_ps( s01, s23, _MM_SHUFFLE( 2, 0, 2, 0 ) );
const __m128 v1 = _mm_shuffle_ps( s01, s23, _MM_SHUFFLE( 3, 1, 3, 1 ) );
// Finally, linear interpolation between these vectors.
const __m128 diff = _mm_sub_ps( v1, v0 );
return _mm_add_ps( v0, _mm_mul_ps( frac, diff ) );
}

程序集看起来不错。现代编译器甚至可以在可用时自动使用 FMA。 (默认情况下，GCC 与-ffp-contract=fast一起执行跨 C 语句的收缩，而不仅仅是在一个表达式中。

---

刚刚看到更新。考虑将目标切换到 SSE 4.1。Steam硬件调查显示，市场渗透率为98.76%。如果您仍然支持像奔腾 4 这样的史前 CPU，那么_mm_floor_ps的解决方法是在 DirectXMath 中，而不是_mm_extract_epi32您可以使用_mm_srli_si128+_mm_cvtsi128_si32。

即使你需要支持像SSE3这样的旧基线，-mtune=generic甚至-mtune=haswell可能是一个好主意，以及-march=nocona，仍然做出内联和其他代码生成选择，这些选择适用于一系列CPU，而不仅仅是奔腾4。