代码之家  ›  专栏  ›  技术社区  ›  markzzz

在SSE2中使用%

  •  3
  • markzzz  · 技术社区  · 7 年前

    下面是我要转换为SSE2的代码:

    double *pA = a;
    double *pB = b[voiceIndex];
    double *pC = c[voiceIndex];
    double *left = audioLeft;
    double *right = audioRight;
    double phase = 0.0;
    double bp0 = mNoteFrequency * mHostPitch;
    
    for (int sampleIndex = 0; sampleIndex < blockSize; sampleIndex++) {
        // some other code (that will use phase)
    
        phase += std::clamp(mRadiansPerSample * (bp0 * pB[sampleIndex] + pC[sampleIndex]), 0.0, PI);
    
        while (phase >= TWOPI) { phase -= TWOPI; }
    }
    

    我的成就如下:

    double *pA = a;
    double *pB = b[voiceIndex];
    double *pC = c[voiceIndex];
    double *left = audioLeft;
    double *right = audioRight;
    double phase = 0.0;
    double bp0 = mNoteFrequency * mHostPitch;
    
    __m128d v_boundLower = _mm_set1_pd(0.0);
    __m128d v_boundUpper = _mm_set1_pd(PI);
    __m128d v_bp0 = _mm_set1_pd(bp0);
    __m128d v_radiansPerSample = _mm_set1_pd(mRadiansPerSample);
    
    __m128d v_phase = _mm_set1_pd(phase);
    __m128d v_pB = _mm_load_pd(pB);
    __m128d v_pC = _mm_load_pd(pC);
    __m128d v_result = _mm_mul_pd(v_bp0, v_pB);
    v_result = _mm_add_pd(v_result, v_pC);
    v_result = _mm_mul_pd(v_result, v_radiansPerSample);
    v_result = _mm_max_pd(v_result, v_boundLower);
    v_result = _mm_min_pd(v_result, v_boundUpper);
    
    for (int sampleIndex = 0; sampleIndex < roundintup8(blockSize); sampleIndex += 8, pB += 8, pC += 8) {
        // some other code (that will use v_phase)
    
        v_phase = _mm_add_pd(v_phase, v_result);
    
        v_pB = _mm_load_pd(pB + 2);
        v_pC = _mm_load_pd(pC + 2);
        v_result = _mm_mul_pd(v_bp0, v_pB);
        v_result = _mm_add_pd(v_result, v_pC);
        v_result = _mm_mul_pd(v_result, v_radiansPerSample);
        v_result = _mm_max_pd(v_result, v_boundLower);
        v_result = _mm_min_pd(v_result, v_boundUpper);
        v_phase = _mm_add_pd(v_phase, v_result);
    
        v_pB = _mm_load_pd(pB + 4);
        v_pC = _mm_load_pd(pC + 4);
        v_result = _mm_mul_pd(v_bp0, v_pB);
        v_result = _mm_add_pd(v_result, v_pC);
        v_result = _mm_mul_pd(v_result, v_radiansPerSample);
        v_result = _mm_max_pd(v_result, v_boundLower);
        v_result = _mm_min_pd(v_result, v_boundUpper);
        v_phase = _mm_add_pd(v_phase, v_result);
    
        v_pB = _mm_load_pd(pB + 6);
        v_pC = _mm_load_pd(pC + 6);
        v_result = _mm_mul_pd(v_bp0, v_pB);
        v_result = _mm_add_pd(v_result, v_pC);
        v_result = _mm_mul_pd(v_result, v_radiansPerSample);
        v_result = _mm_max_pd(v_result, v_boundLower);
        v_result = _mm_min_pd(v_result, v_boundUpper);
        v_phase = _mm_add_pd(v_phase, v_result);
    
        v_pB = _mm_load_pd(pB + 8);
        v_pC = _mm_load_pd(pC + 8);
        v_result = _mm_mul_pd(v_bp0, v_pB);
        v_result = _mm_add_pd(v_result, v_pC);
        v_result = _mm_mul_pd(v_result, v_radiansPerSample);
        v_result = _mm_max_pd(v_result, v_boundLower);
        v_result = _mm_min_pd(v_result, v_boundUpper);
    
        // ... fmod?
    }
    

    但我真的不知道怎么换 while (phase >= TWOPI) { phase -= TWOPI; } (基本上是经典 fmod 在C++中。

    有什么奇特的内在因素吗?在这上面找不到任何 list . 除法+某种火箭钻头移动?

    1 回复  |  直到 7 年前
        1
  •  4
  •   llllllllll    7 年前

    正如评论所说,在这个例子中,你可以把它变成一个带比较的蒙面减法。+ andpd . 只要从返回到所需的范围中减去一个以上,这就有效。

    喜欢

    const __m128d v2pi = _mm_set1_pd(TWOPI);
    
    
    __m128d needs_range_reduction = _mm_cmpge_pd(vphase, v2pi);
    __m128d offset = _mm_and_pd(needs_range_reduction, v2pi);  // 0.0 or 2*Pi
    vphase = _mm_sub_pd(vphase, offset);
    

    实现一个实际的(缓慢的) fmod 在不太担心最后几点意义的情况下,你会这样做的。 integer_quotient = floor(x/y) (或者) rint(x/y) 或 ceil ) x - y * integer_quotient . floor / rint / 塞尔 SSE4.1便宜 _mm_round_pd 或 _mm_floor_pd() . 这将给出余数,它可以是负数,就像整数除法一样。

    我相信有一些数值技术可以更好地避免在灾难性的取消之前从两个附近的数字中减去舍入误差。如果你关心精度,去检查一下。(使用) double 向量当你不太关心精确性的时候是有点愚蠢的;不妨使用 float 得到每个向量两倍的工作量)。如果输入比模大得多,则会不可避免地损失精度,而最小化临时舍入误差可能非常重要。但否则,除非您关心结果中非常接近零的相对误差,否则精度只会是一个问题。 x 几乎是 y . (接近零的结果,只剩下有效位的底端几位用于精确。)

    如果没有SSE4.1,有一些技巧,比如加上然后减去一个足够大的数字。转换为整数和后整数对于 pd 因为压缩的转换指令也会解码为一些无序的UOP。更不用说32位整数不能覆盖 双重的 但是,如果你的输入量那么大,你就无法获得精确的减程。

    如果你有 FMA ,可以避免 y * integer_quotient 乘法和Sub的一部分。 _mm_fmsub_pd .

    推荐文章