一些mandelbrot绘制例程从c到sse2

我想重写这样简单的例程到SSE2代码,(最好是在nasm中)我并不完全确定如何做,两件事情不清楚(如何表达计算(内循环和那些来自外循环)以及如何调用c代码函数“SetPixelInDibInt(i,j,palette [n]);”来自staticaly链接的asm代码

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 4.0 ) break; n++; } SetPixelInDibInt(i ,j, palette[n]); } } } 

可能有人帮助,我不希望看到其他代码实现,但只是上面那些的鼻子翻译 – 在我的情况下,这将是最有帮助的 – 有人可以帮忙吗?

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

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

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

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

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

可以使用打包比较来检查循环条件,以生成零和1的掩码,以及该寄存器的(V)PTEST以查看它是否全为零。 (编辑: movmskps然后test+jcc是更少的movmskps ,但可能更高的延迟。)然后显然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上使用PTEST ,如果FP值大于4.0,则会选择一些位错误的AND-mask来选择将被设置的位。 就像指数中的一些比特一样? 也许值得考虑。 我发现了一个关于它的论坛post ,但没试过。

嗯,哦,废话,这不记录循环条件失败,为每个矢量元素分别,为了着色Mandelbrot集外的点。 也许测试任何符合条件的元素(而不是全部),记录结果,然后将该元素(和该元素的c )设置为0.0,这样它就不会再次触发退出条件。 或者可能将像素调度到矢量元素中是最重要的。 这个代码在超线程CPU上可能做得相当好,因为会有很多分支错误预测,每个元素分别触发早期状态。

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

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

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

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

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

他们只是让FP结果不断累积到早期失败状态的元素。 他们只是停止递增该元素的计数器。 希望大多数系统默认将控制字设置为零非正规,如果非正规仍然需要额外的周期。

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

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