把复数数组的乘积向量化

Question

我正在尝试编写快速/最优的代码，将复数数组的乘积矢量化。在简单的C中，这将是：

#include <complex.h>
complex float f(complex float x[], int n ) {
  complex float p = 1.0;
  for (int i = 0; i < n; i++)
    p *= x[i];
  return p;
}

然而，gcc不能矢量化这一点，我的目标CPU是AMD FX-8350，它支持AVX。为了加快速度，我试过：

typedef float v4sf __attribute__ ((vector_size (16)));
typedef union {
  v4sf v;
  float e[4];
} float4;
typedef struct {
  float4 x;
  float4 y;
} complex4;
static complex4 complex4_mul(complex4 a, complex4 b) {
  return (complex4){a.x.v*b.x.v -a.y.v*b.y.v, a.y.v*b.x.v + a.x.v*b.y.v};
}
complex4 f4(complex4 x[], int n) {
  v4sf one = {1,1,1,1};
  complex4 p = {one,one};
  for (int i = 0; i < n; i++) p = complex4_mul(p, x[i]);
  return p;
}

对于AVX，这段代码可以改进吗?它是否尽可能快？

Answer 1

我想知道编译器是否能够更容易地将其矢量化，而不需要union内的struct和e成员的潜在类型双关语(这似乎没有被使用)。

这个怎么样？

typedef float v4sf __attribute__ ((vector_size (16)));
typedef struct {
  v4sf x;
  v4sf y;
} complex4;

static inline v4sf complex4_mul_r(v4sf a_r, v4sf a_i, v4sf b_r, v4sf b_i) {
  return a_r*b_r -a_i*b_i;
}
static inline v4sf complex4_mul_i(v4sf a_r, v4sf a_i, v4sf b_r, v4sf b_i) {
  return a_r*b_i + a_i*b_r;
}

complex4 f4(v4sf x_r[], v4sf x_i[], int n) {
  v4sf one = {1,1,1,1};
  v4sf p_r = one;
  v4sf p_i = one;
  v4sf p_r_temp;
  for (int i = 0; i < n; i++)
  {
     p_r_temp = complex4_mul_r(p_r, p_i, x_r[i], x_i[i]);
     p_i = complex4_mul_i(p_r, p_i, x_r[i], x_i[i]);
     p_r = p_r_temp;
  }
  return (complex4){p_r, p_i};
}

看一下https://godbolt.org/上的程序集，它似乎得到了完全的矢量化。我无法让“上帝”共享链接开始工作。

我在考虑是否有可能将实部和虚部放在同一个向量中，并在必要时使用__builtin_shuffle()对它们进行排列。不能完全解决这个问题。