mandelbrot的simd实现方式
我用普通方式实现了下面的函数,计算Mandelbrot集合的数量积
// --------------------------------------------------
size_t mandelbrot_scalar(float a, float b, size_t max_iter)
// --------------------------------------------------
{
size_t num_iter = 0;
float x,t,t2,y;
x = 0.0f;
y = 0.0f;
t = 0.0f;
t2= 0.0f;
//iteration 0
y= 2*x*y+b;
x=t-t2+a;
t=x*x;
t2=y*y;
//calcul
while((num_iter < max_iter) && (t <= 4.0) ) {
t=x*x;
t2=y*y;
y= 2*x*y+b;
x=t-t2+a;
t = t+t2;
num_iter++;
};
//cas non iter max
if(num_iter < max_iter)
num_iter--;
return num_iter;
}
现在我想用simd的方式实现以上代码,我的想法是用simd指令集中对应的函数来替换以上函数中的加减乘除以及比较等运算符号
// --------------------------------------------------------------
vuint32 mandelbrot_simd(vfloat32 a, vfloat32 b, size_t max_iter)
// --------------------------------------------------------------
{
vuint32 num_iter = _mm_set1_epi32(0);
vfloat32 x = _mm_setzero_ps();
vfloat32 y = _mm_setzero_ps();
vfloat32 t = _mm_setzero_ps();
vfloat32 t2 = _mm_setzero_ps();
vfloat32 cons0 = _mm_set1_ps(0.f);
vuint32 cons1 = _mm_set1_epi32(1.f);
vfloat32 cons2 = _mm_set1_ps(2.f);
vfloat32 cons4 = _mm_set1_ps(4.f);
y=_mm_add_ps(_mm_mul_ps(_mm_mul_ps(cons2,x),y),b);
x=_mm_add_ps(_mm_sub_ps(t,t2),a);
t=_mm_mul_ps(x,x);
t2=_mm_mul_ps(y,y);
while(!(_mm_movemask_ps((_mm_cmplt_ps(_mm_cvtepi32_ps(num_iter),_mm_cvtepi32_ps(_mm_set1_epi32(max_iter)))))+1)&&!(_mm_movemask_ps((_mm_cmple_ps(t,cons4)))+1)) {
t=_mm_mul_ps(x,x);
t2=_mm_mul_ps(y,y);
y=_mm_add_ps(_mm_mul_ps(_mm_mul_ps(cons2,x),y),b);
x=_mm_add_ps(_mm_sub_ps(t,t2),a);
t=_mm_add_ps(t,t2);
num_iter=_mm_add_epi32(num_iter,cons1);
};
//cas non iter max
if(!(_mm_movemask_ps((_mm_cmplt_ps(_mm_cvtepi32_ps(num_iter),_mm_cvtepi32_ps(_mm_set1_epi32(max_iter))))))+1)
num_iter=_mm_sub_epi32(num_iter,cons1);
return num_iter;
}
但是严格替换下来最终结果不对,不知道哪里出了问题,有大神指导下吗