2

I want to rewrite such simple routine to SSE2 code, (preferably in nasm) and I am not totally sure how to do it, two things not clear (how to express calculations (inner loop and those from outer loop too) and how to call c code function "SetPixelInDibInt(i ,j, palette[n]);" from under staticaly linked asm code

    void DrawMandelbrotD(double ox, double oy, double lx, int N_ITER)
    {
     double ly = lx * double(CLIENT_Y)/double(CLIENT_X);
     double dx = lx / CLIENT_X;
     double dy = ly / CLIENT_Y;
     double ax = ox - lx * 0.5 + dx * 0.5;
     double ay = oy - ly * 0.5 + dy * 0.5;
    static  double re, im, re_n, im_n, c_re, c_im, rere, imim, int n;

    for(int j=0; j<CLIENT_Y; j+=1)
    {
     for(int i=0; i<CLIENT_X; i+=1)
     {
      c_re = ax + i * dx;
      c_im = ay + j * dy;
      re = c_re;
      im = c_im;
      rere=re*re;
      imim=im*im;
      n=1;

      for(int k=0;k<N_ITER;k++)
      {
        im =  (re+re)*im    + c_im;
        re =   rere - imim  + c_re;
        rere=re*re;
        imim=im*im;
        if ( (rere + imim) > 4.0 ) break;
        n++;
       }
        SetPixelInDibInt(i ,j, palette[n]);
      }
     }
    }

could someone help, I would like not to see other code implementations but just nasm-sse translation of those above - it would be most helpfull in my case - could someone help with that?

4

1 回答 1

2

英特尔有一个完整的实现作为 AVX 示例。见下文。

使 Mandelbrot 棘手的是集合中每个点(即像素)的早期输出条件不同。您可以保持一对或四倍像素的迭代,直到两者的大小都超过 2.0(或者您达到最大迭代次数)。否则需要跟踪哪个像素的点在哪个向量元素中。

无论如何,一次对 2 个(或 4 个使用 AVX)双精度向量进行操作的简单实现将使其吞吐量受到依赖链延迟的限制。您需要并行执行多个依赖链,以保持 Haswell 的两个 FMA 单元的供给。因此,您将复制变量,并在内部循环内对外部循环的两次迭代进行交错操作。

跟踪正在计算的像素会有点棘手。我认为为一行像素使用一组寄存器,为另一行使用另一组寄存器可能需要更少的开销。(所以你总是可以向右移动 4 个像素,而不是检查其他 dep 链是否已经在处理该向量。)

我怀疑每 4 次左右迭代只检查一次循环退出条件可能是一个胜利。基于压缩向量比较使代码分支,比标量情况稍微贵一些。所需的额外 FP 添加也很昂贵。(Haswell 每个周期可以执行两个 FMA,(延迟 = 5)。唯一的 FP add 单元与 FMA 单元之一是同一端口。两个 FP mul 单元位于可以运行 FMA 的相同端口上。)

可以使用打包比较来检查循环条件以生成零和一的掩码,以及(V)PTEST该寄存器的一个自身以查看它是否全为零。(编辑:movmskps然后test+jcc是更少的微指令,但可能会更高的延迟。)然后显然jejne适当,这取决于您是否进行了 FP 比较,在您应该退出时留下零,或者在您不应该退出时留下零。NAN 不应该是可能的,但是没有理由不选择您的比较操作,这样 NAN 将导致退出条件为真。

const __mm256d const_four = _mm256_set1_pd(4.0);  // outside the loop

__m256i cmp_result = _mm256_cmp_pd(mag_squared, const_four, _CMP_LE_OQ);  // vcmppd.  result is non-zero if at least one element < 4.0
if (_mm256_testz_si256(cmp_result, cmp_result))
    break;

可能有某种方法可以PTEST直接在打包双精度上使用,并带有一些位黑客 AND 掩码,它将挑选出将在 FP 值 > 4.0 时设置的位。就像指数中的一些位?也许值得考虑。我找到了一个关于它的论坛帖子,但没有尝试过。

嗯,哦,废话,这不会记录循环条件失败的时间,分别针对每个向量元素,目的是为 Mandelbrot 集之外的点着色。也许测试任何符合条件的元素(而不是全部),记录结果,然后将该元素(以及c该元素)设置为 0.0,这样它就不会再次触发退出条件。或者,也许将像素安排到矢量元素中才是最终的方法。这段代码在超线程 CPU 上可能表现得相当好,因为每个元素都会单独触发提前退出条件,因此会出现很多分支错误预测。

这可能会浪费您的大量吞吐量,并且考虑到每个周期 4 uops 是可行的,但其中只有 2 个可以是 FP mul/add/FMA,因此有大量整数代码将点安排到向量元素中的空间。(在 Sandybridge/Ivybrideg 上,没有 FMA,FP 吞吐量较低。但只有 3 个端口可以处理整数操作,其中 2 个是 FP mul 和 FP add 单元的端口。)

由于您不必读取任何源数据,因此每个 dep 链只有 1 个内存访问流,它是一个写入流。(而且它的带宽很低,因为大多数点在您准备好写入单个像素值之前需要进行大量迭代。)因此硬件预取流的数量并不是并行运行的 dep 链数量的限制因素. 缓存未命中延迟应被写入缓冲区隐藏。

如果有人对此仍然感兴趣,我可以编写一些代码(只需发表评论)。不过,我在高级设计阶段停了下来,因为这是一个老问题。

===============

我还发现英特尔已经使用 Mandelbrot 集作为其AVX 教程之一的示例。他们对循环条件使用 mask-off-vector-elements 方法。(使用vcmppsto AND 直接生成的掩码)。他们的结果表明,AVX(单精度)比标量浮点数提高了 7 倍,因此显然相邻像素在非常不同的迭代次数下达到提前输出条件并不常见。(至少对于他们测试的缩放/平移。)

他们只是让 FP 结果不断累积,以获取未能满足早期退出条件的元素。他们只是停止增加该元素的计数器。希望大多数系统默认将控制字设置为零,如果非正规仍然需要额外的周期。

但是,他们的代码在某一方面很愚蠢:他们使用浮点向量跟踪每个向量元素的迭代计数,然后在使用前将其转换为 int。使用压缩整数会更快,并且不占用 FP 执行单元。哦,我知道他们为什么这样做:AVX(没有 AVX2)不支持 256 位整数向量操作。他们本可以使用打包的 16 位 int 循环计数器,但这可能会溢出。(而且他们必须将掩码从 256b 压缩到 128b)。

他们还测试所有 > 4.0 的元素movmskps,然后测试,而不是使用ptest. 我猜test / jcccan 宏融合,并在与 FP 矢量操作不同的执行单元上运行,所以它可能不会更慢。哦,当然 AVX(没有 AVX2)没有 256bit PTEST。此外,PTEST是 2 微指令,所以实际上movmskps+比 .test / jcc微指令少ptest + jcc。(PTEST在 SnB 上是 1 个融合域 uop,但对于执行端口仍有 2 个未融合的 uop。在 IvB/HSW 上,即使在融合域中也是 2 个 uops。)所以它看起来movmskps是最佳方式,除非您可以利用按位 AND 这是 的一部分PTEST,或者需要测试的不仅仅是每个元素的高位。如果一个分支是不可预测的,ptest可能会降低延迟,因此通过更快地捕获错误预测一个周期是值得的。

于 2015-06-25T21:23:12.567 回答