使用 SIMD 将累积（单个）值打包成两个值管理清理代码循环的方法是啥？

Posted 2023-02-16

【中文标题】使用 SIMD 将累积（单个）值打包成两个值管理清理代码循环的方法是啥？

【英文标题】：What the method to manage Cleanup Code loop for a cumulative (single) value packed into two values using SIMD?使用 SIMD 将累积（单个）值打包成两个值管理清理代码循环的方法是什么？

【发布时间】：2019-01-10 09:55:56

【问题描述】：

假设我管理一个名为v_phase 的__m128d 变量，其计算方式为

index 0 : load prev phase
index 1 : phase += newValue
index 2 : phase += newValue
index 3 : phase += newValue
index 4 : phase += newValue
...

这是基本代码：

__m128d v_phase;

// load prev cumulated mPhase to v_phase (as mPhase, mPhase + nextValue)

for (int sampleIndex = 0; sampleIndex < blockSize; sampleIndex += 2, pValue += 2) 
    // function with phase

    // update pValue increment (its not linear)

    // phase increment: v_phase += newValue


// cleanup code
if (blockSize % 2 == 0) 
    mPhase = v_phase.m128d_f64[0];

事实是：如果blockSize 是偶数，它工作正常：它将在最后一次循环迭代中求和两个相位值，并取v_phase.m128d_f64[0]（即第一个新增两个）。

但是如果blockSize 是奇数呢？我只需要最后一次迭代的v_phase.m128d_f64[1]不加两个相位值。

我可以使用sampleIndex < blockSize - 1，但这会在// cleanup code 中移动逻辑// function with phase（我不太喜欢它）。

在循环中放置一个 if 是我会避免的（branc 预测；因为我使用的是 SIMD，所以我正在优化代码，这会变慢）。

有什么建议吗？

这是一个更“完整”的例子：

double phase = mPhase;

__m128d v_pB = _mm_setr_pd(0.0, pB[0]);
v_pB = _mm_mul_pd(v_pB, v_radiansPerSampleBp0);
__m128d v_pC = _mm_setr_pd(0.0, pC[0]);
v_pC = _mm_mul_pd(v_pC, v_radiansPerSample);

__m128d v_pB_prev = _mm_setr_pd(0.0, 0.0);
v_pB_prev = _mm_mul_pd(v_pB_prev, v_radiansPerSampleBp0);
__m128d v_pC_prev = _mm_setr_pd(0.0, 0.0);
v_pC_prev = _mm_mul_pd(v_pC_prev, v_radiansPerSample);

__m128d v_phaseAcc1;
__m128d v_phaseAcc2;
__m128d v_phase = _mm_set1_pd(phase);

// phase
v_phaseAcc1 = _mm_add_pd(v_pB, v_pC);
v_phaseAcc1 = _mm_max_pd(v_phaseAcc1, v_boundLower);
v_phaseAcc1 = _mm_min_pd(v_phaseAcc1, v_boundUpper);
v_phaseAcc2 = _mm_add_pd(v_pB_prev, v_pC_prev);
v_phaseAcc2 = _mm_max_pd(v_phaseAcc2, v_boundLower);
v_phaseAcc2 = _mm_min_pd(v_phaseAcc2, v_boundUpper);
v_phase = _mm_add_pd(v_phase, v_phaseAcc1);
v_phase = _mm_add_pd(v_phase, v_phaseAcc2);

for (int sampleIndex = 0; sampleIndex < blockSize; sampleIndex += 2, pB += 2, pC += 2) 
    // code that will use v_phase

    // phase increment
    v_pB = _mm_loadu_pd(pB + 1);
    v_pB = _mm_mul_pd(v_pB, v_radiansPerSampleBp0);
    v_pC = _mm_loadu_pd(pC + 1);
    v_pC = _mm_mul_pd(v_pC, v_radiansPerSample);

    v_pB_prev = _mm_load_pd(pB);
    v_pB_prev = _mm_mul_pd(v_pB_prev, v_radiansPerSampleBp0);
    v_pC_prev = _mm_load_pd(pC);
    v_pC_prev = _mm_mul_pd(v_pC_prev, v_radiansPerSample);

    v_phaseAcc1 = _mm_add_pd(v_pB, v_pC);
    v_phaseAcc1 = _mm_max_pd(v_phaseAcc1, v_boundLower);
    v_phaseAcc1 = _mm_min_pd(v_phaseAcc1, v_boundUpper);
    v_phaseAcc2 = _mm_add_pd(v_pB_prev, v_pC_prev);
    v_phaseAcc2 = _mm_max_pd(v_phaseAcc2, v_boundLower);
    v_phaseAcc2 = _mm_min_pd(v_phaseAcc2, v_boundUpper);
    v_phase = _mm_add_pd(v_phase, v_phaseAcc1);
    v_phase = _mm_add_pd(v_phase, v_phaseAcc2);


// cleanup code
if (blockSize % 2 == 0) 
    mPhase = v_phase.m128d_f64[0];

else 
    ??? if odd?

【问题讨论】：

也许只有我一个人，但问题似乎不清楚。例如，“pValue increment”是什么意思？您可能想添加一个minimal reproducible example。

我可以给出的一般建议是 (1) 重新排列循环部分和 (2) 将阶段初始化为较低的值，对应于迭代 -1。但是如果不查看所有详细信息，就无法确定它是否有用。

"将阶段初始化为较低的值，对应于迭代-1"我不能。总有一个起点，它的“0”、-1 不存在

@anatolyg：我已经添加了更多代码，现在应该清楚了！

【参考方案1】：

除了最后一个之外，您还可以从循环中输出 previous v_phase。也就是说，在更新你的v_phase之前，先保存上一个：

__m128d prev_v_phase;
for (...) 
    ...
    prev_v_phase = v_phase;
    v_phase = _mm_add_pd(v_phase, v_phaseAcc1);
    v_phase = _mm_add_pd(v_phase, v_phaseAcc2);


// cleanup code
if (blockSize % 2 == 0) 
    mPhase = v_phase.m128d_f64[0];

else 
    mPhase = prev_v_phase.m128d_f64[1];

如果循环根本不执行任何迭代，这将失败（然后prev_v_phase 将未初始化），但这是性能不重要的情况，很容易处理。

【讨论】：

它似乎有效。对于一个常见的价值增量来说，可能有点“太复杂”了？ :)