SIMD寄存器的数学函数 [英] Mathematical functions for SIMD registers
问题描述
根据 https://sourceware.org/glibc/wiki/libmvec GCC拥有数学函数的向量实现.可以在编译器中使用它们进行优化，在本示例中可以看到: https://godbolt.org/g/IcxtVi ，编译器会使用一些错误的正弦函数，并且一次操作4个双精度数
According to https://sourceware.org/glibc/wiki/libmvec GCC has vector implementation of math functions. They can be used by compiler for optimizations, it can be seen in this example: https://godbolt.org/g/IcxtVi, compiler uses some mangled sine function and operates on 4 doubles at a time
我知道有一些SIMD数学库可以在需要数学函数的情况下使用，但我仍然感兴趣的是，有一种方法可以通过某种方式手动调用__m256d变量中GCC中已经存在的矢量化数学函数.固有的还是其他方式?
I know that there are SIMD math libraries that can be used if I need math functions, but I am still interested is there a way to manually call vectorized math functions that already exist in GCC on __m256d variable with some kind of intrinsic or in any other way?
推荐答案
Libmvec是一个x86_64 glibc库，具有用于以下目的的SSE4，AVX，AVX2和AVX-512矢量化函数
cos，exp，log，sin，pow和sincos，具有单精度和双精度.
这些函数的精度是最大4ulp的相对误差.
Libmvec is a x86_64 glibc library with SSE4, AVX, AVX2, and AVX-512 vectorized functions for
cos, exp, log, sin, pow, and sincos, in single precision and double precision.
The accuracy of these functions is 4-ulp maximum relative error.
通常，gcc会插入对Libmvec函数的调用，例如_ZGVdN4v_cos，
当它自动矢量化标量代码时，例如使用-ffast-math或
#pragma omp simd.但是，这些向量函数也适用于
手动向量化C代码.
不幸的是， Libmvec文档
没有给出很多示例来说明如何将这些功能用于
手动向量化C代码.
Usually gcc inserts calls to Libmvec functions, such as _ZGVdN4v_cos,
while it is auto-vectorizing scalar code, for example with -ffast-math or
#pragma omp simd. However, these vector functions are also suitable for
manually vectorizing C code.
Unfortunately the Libmvec document
does not give many examples that explain how to use these functions for
manually vectorizing C code.
要使用Libmvec函数，必须包含函数声明
并用-lm(与数学库链接)选项链接代码.
链接器不需要-lmvec，链接标准数学库-lm是
对于非静态版本而言已足够，请参见此处.
-ffast-math编译器选项不是必需的.
To use the Libmvec functions one has to include the function declarations
and link the code with the -lm (link with math library) option.
The linker doesn't need -lmvec, linking the standard math library -lm is
sufficient for not static builds, see here.
The -ffast-math compiler option is not necesary.
以下示例显示如何调用AVX2矢量化函数. 此示例可在标准Ubuntu 18.04系统(glibc版本2.27)上成功运行.
The following example shows how to call the AVX2 vectorized functions. This example runs successful on a standard Ubuntu 18.04 system (glibc version 2.27).
/* gcc -Wall -O3 -m64 -march=skylake libmvec_ex.c -lm */
#include <immintrin.h>
#include <stdio.h>
#include <stdint.h>
#define _GNU_SOURCE /* These two lines are optional. They are only needed to use the scalar */
#include <math.h> /* functions sin, exp, sincos,... */
__m256d _ZGVdN4v_cos(__m256d x);
__m256d _ZGVdN4v_exp(__m256d x);
__m256d _ZGVdN4v_log(__m256d x);
__m256d _ZGVdN4v_sin(__m256d x);
__m256d _ZGVdN4vv_pow(__m256d x, __m256d y);
void _ZGVdN4vvv_sincos(__m256d x, __m256i ptrs, __m256i ptrc);
__m256 _ZGVdN8v_cosf(__m256 x);
__m256 _ZGVdN8v_expf(__m256 x);
__m256 _ZGVdN8v_logf(__m256 x);
__m256 _ZGVdN8v_sinf(__m256 x);
__m256 _ZGVdN8vv_powf(__m256 x, __m256 y);
void _ZGVdN8vvv_sincosf(__m256 x, __m256i ptrs_lo, __m256i ptrs_hi, __m256i ptrc_lo, __m256i ptrc_hi);
int mm256_print_pd(__m256d x);
int mm256_print_ps(__m256 x);
int main()
{
__m256d x, y, z;
__m256 xf, yf, zf;
double z_s[4], z_c[4]; /* sincos output */
float zf_s[8], zf_c[8]; /* sincosf output */
__m256i ptrs, ptrc; /* Pointers to the elements of z_s and z_c */
__m256i ptrs_lo, ptrs_hi, ptrc_lo, ptrc_hi;
x = _mm256_set_pd (0.04, 0.03, 0.02, 0.01);
y = _mm256_set_pd (2.0, 1.5, 1.0, 0.5);
xf = _mm256_set_ps (0.08f, 0.07f, 0.06f, 0.05f, 0.04f, 0.03f, 0.02f, 0.01f);
yf = _mm256_set_ps (4.0f, 3.5f, 3.0f, 2.5f, 2.0f, 1.5f, 1.0f, 0.5f);
printf("AVX2 Double precision examples\n");
printf("x "); mm256_print_pd(x);
printf("y "); mm256_print_pd(y);
z =_ZGVdN4v_cos(x); printf("cos(x) "); mm256_print_pd(z);
z =_ZGVdN4v_exp(x); printf("exp(x) "); mm256_print_pd(z);
z =_ZGVdN4v_log(x); printf("log(x) "); mm256_print_pd(z);
z =_ZGVdN4v_sin(x); printf("sin(x) "); mm256_print_pd(z);
z =_ZGVdN4vv_pow(x, y); printf("pow(x,y) "); mm256_print_pd(z);
ptrs = _mm256_set_epi64x((uint64_t)&z_s[3],(uint64_t)&z_s[2],(uint64_t)&z_s[1],(uint64_t)&z_s[0]);
ptrc = _mm256_set_epi64x((uint64_t)&z_c[3],(uint64_t)&z_c[2],(uint64_t)&z_c[1],(uint64_t)&z_c[0]);
/* Alternative: ptrs = _mm256_add_epi64(_mm256_set1_epi64x((uint64_t)&z_s[0]),_mm256_set_epi64x(24,16,8,0)); */
/* This might be more efficient if the destination addresses are contiguous in memory. */
_ZGVdN4vvv_sincos(x, ptrs, ptrc); /* The results of _ZGVdN4vvv_sincos are scattered into the adresses in ptrs and ptrc */
printf("sincos cos(x) %12.8f %12.8f %12.8f %12.8f \n", z_c[3], z_c[2], z_c[1], z_c[0]);
printf("sincos sin(x) %12.8f %12.8f %12.8f %12.8f\n\n", z_s[3], z_s[2], z_s[1], z_s[0]);
printf("AVX2 Single precision examples\n");
printf("x "); mm256_print_ps(xf);
printf("y "); mm256_print_ps(yf);
zf =_ZGVdN8v_cosf(xf); printf("cosf(x) "); mm256_print_ps(zf);
zf =_ZGVdN8v_expf(xf); printf("expf(x) "); mm256_print_ps(zf);
zf =_ZGVdN8v_logf(xf); printf("logf(x) "); mm256_print_ps(zf);
zf =_ZGVdN8v_sinf(xf); printf("sinf(x) "); mm256_print_ps(zf);
zf =_ZGVdN8vv_powf(xf, yf);printf("powf(x,y) "); mm256_print_ps(zf);
ptrs_lo = _mm256_set_epi64x((uint64_t)&zf_s[3],(uint64_t)&zf_s[2],(uint64_t)&zf_s[1],(uint64_t)&zf_s[0]);
ptrs_hi = _mm256_set_epi64x((uint64_t)&zf_s[7],(uint64_t)&zf_s[6],(uint64_t)&zf_s[5],(uint64_t)&zf_s[4]);
ptrc_lo = _mm256_set_epi64x((uint64_t)&zf_c[3],(uint64_t)&zf_c[2],(uint64_t)&zf_c[1],(uint64_t)&zf_c[0]);
ptrc_hi = _mm256_set_epi64x((uint64_t)&zf_c[7],(uint64_t)&zf_c[6],(uint64_t)&zf_c[5],(uint64_t)&zf_c[4]);
_ZGVdN8vvv_sincosf(xf, ptrs_lo, ptrs_hi, ptrc_lo, ptrc_hi); /* The results of _ZGVdN8vvv_sincosf are scattered to the adresses in ptrs_lo, ptrs_hi, ptrc_lo, and ptrc_hi */
printf("sincosf cos(x)%12.8f %12.8f %12.8f %12.8f %12.8f %12.8f %12.8f %12.8f \n",
zf_c[7], zf_c[6], zf_c[5], zf_c[4], zf_c[3], zf_c[2], zf_c[1], zf_c[0]);
printf("sincosf sin(x)%12.8f %12.8f %12.8f %12.8f %12.8f %12.8f %12.8f %12.8f \n",
zf_s[7], zf_s[6], zf_s[5], zf_s[4], zf_s[3], zf_s[2], zf_s[1], zf_s[0]);
return 0;
}
__attribute__ ((noinline)) int mm256_print_pd(__m256d x){
double vec_x[4];
_mm256_storeu_pd(vec_x,x);
printf("%12.8f %12.8f %12.8f %12.8f \n", vec_x[3], vec_x[2], vec_x[1], vec_x[0]);
return 0;
}
__attribute__ ((noinline)) int mm256_print_ps(__m256 x){
float vec_x[8];
_mm256_storeu_ps(vec_x,x);
printf("%12.8f %12.8f %12.8f %12.8f %12.8f %12.8f %12.8f %12.8f \n", vec_x[7], vec_x[6], vec_x[5], vec_x[4],
vec_x[3], vec_x[2], vec_x[1], vec_x[0]);
return 0;
}
/*
输出符合预期:
$ ./a.out
AVX2 Double precision examples
x 0.04000000 0.03000000 0.02000000 0.01000000
y 2.00000000 1.50000000 1.00000000 0.50000000
cos(x) 0.99920011 0.99955003 0.99980001 0.99995000
exp(x) 1.04081077 1.03045453 1.02020134 1.01005017
log(x) -3.21887582 -3.50655790 -3.91202301 -4.60517019
sin(x) 0.03998933 0.02999550 0.01999867 0.00999983
pow(x,y) 0.00160000 0.00519615 0.02000000 0.10000000
sincos cos(x) 0.99920011 0.99955003 0.99980001 0.99995000
sincos sin(x) 0.03998933 0.02999550 0.01999867 0.00999983
AVX2 Single precision examples
x 0.08000000 0.07000000 0.06000000 0.05000000 0.04000000 0.03000000 0.02000000 0.01000000
y 4.00000000 3.50000000 3.00000000 2.50000000 2.00000000 1.50000000 1.00000000 0.50000000
cosf(x) 0.99680173 0.99755102 0.99820060 0.99875027 0.99920011 0.99955004 0.99979997 0.99995005
expf(x) 1.08328700 1.07250810 1.06183648 1.05127108 1.04081070 1.03045452 1.02020133 1.01005018
logf(x) -2.52572870 -2.65926003 -2.81341076 -2.99573231 -3.21887589 -3.50655794 -3.91202307 -4.60517025
sinf(x) 0.07991469 0.06994285 0.05996400 0.04997917 0.03998933 0.02999550 0.01999867 0.00999983
powf(x,y) 0.00004096 0.00009075 0.00021600 0.00055902 0.00160000 0.00519615 0.02000000 0.10000000
sincosf cos(x) 0.99680173 0.99755102 0.99820060 0.99875027 0.99920011 0.99955004 0.99979997 0.99995005
sincosf sin(x) 0.07991469 0.06994285 0.05996400 0.04997917 0.03998933 0.02999550 0.01999867 0.00999983
请注意，向量化的sincos函数会将结果分散到4(双精度)或
8(浮动)，不一定是连续的内存位置.
Note that the vectorized sincos functions scatters the results to 4 (double) or
8 (float), not necessarily contiguous, memory locations.
SSE4，AVX和AVX2的函数声明，以及Intel SVML的替代方法.
下面的代码块显示了许多不同情况下的函数声明， 以及英特尔SVML库中的相应SVML函数. 不知何故我找不到带有正确参数列表的这些声明， 在互联网上.
The code block below shows the function declaration for many different cases, and the corresponding SVML function from Intel's SVML library. Somehow I couldn't find these declarations, with the right parameter lists, on the internet.
/* Double precision Libmvec Intel SVML function */
/* SSE4 */
__m128d _ZGVbN2v_cos(__m128d x); /* _mm_cos_pd(x) */
__m128d _ZGVbN2v_exp(__m128d x); /* _mm_exp_pd(x) */
__m128d _ZGVbN2v_log(__m128d x); /* _mm_log_pd(x) */
__m128d _ZGVbN2v_sin(__m128d x); /* _mm_sin_pd(x) */
__m128d _ZGVbN2vv_pow(__m128d x, __m128d y); /* _mm_pow_pd(x, y) */
void _ZGVbN2vvv_sincos(__m128d x, __m128i ptrs, __m128i ptrc); /* _mm_sincos_pd */
/* AVX */
__m256d _ZGVcN4v_cos(__m256d x); /* _mm256_cos_pd(x) */
__m256d _ZGVcN4v_exp(__m256d x); /* _mm256_exp_pd(x) */
__m256d _ZGVcN4v_log(__m256d x); /* _mm256_log_pd(x) */
__m256d _ZGVcN4v_sin(__m256d x); /* _mm256_sin_pd(x) */
__m256d _ZGVcN4vv_pow(__m256d x, __m256d y); /* _mm256_pow_pd(x, y) */
void _ZGVcN4vvv_sincos(__m256d x, __m256i ptrs, __m256i ptrc); /* _mm256_sincos_pd */
/* AVX2 */
__m256d _ZGVdN4v_cos(__m256d x); /* _mm256_cos_pd(x) */
__m256d _ZGVdN4v_exp(__m256d x); /* _mm256_exp_pd(x) */
__m256d _ZGVdN4v_log(__m256d x); /* _mm256_log_pd(x) */
__m256d _ZGVdN4v_sin(__m256d x); /* _mm256_sin_pd(x) */
__m256d _ZGVdN4vv_pow(__m256d x, __m256d y); /* _mm256_pow_pd(x, y) */
void _ZGVdN4vvv_sincos(__m256d x, __m256i ptrs, __m256i ptrc); /* _mm256_sincos_pd */
/* Single precision Libmvec Intel SVML function */
/* SSE4 */
__m128 _ZGVbN4v_cosf(__m128 x); /* _mm_cos_ps(x) */
__m128 _ZGVbN4v_expf(__m128 x); /* _mm_exp_ps(x) */
__m128 _ZGVbN4v_logf(__m128 x); /* _mm_log_ps(x) */
__m128 _ZGVbN4v_sinf(__m128 x); /* _mm_sin_ps(x) */
__m128 _ZGVbN4vv_powf(__m128 x, __m128 y); /* _mm_pow_ps(x, y) */
void _ZGVbN4vvv_sincosf(__m128 x, __m128i ptrs_lo, __m128i ptrs_hi, __m128i ptrc_lo, __m128i ptrc_hi); /* _mm_sincos_ps */
/* AVX */
__m256 _ZGVcN8v_cosf(__m256 x); /* _mm256_cos_ps(x) */
__m256 _ZGVcN8v_expf(__m256 x); /* _mm256_exp_ps(x) */
__m256 _ZGVcN8v_logf(__m256 x); /* _mm256_log_ps(x) */
__m256 _ZGVcN8v_sinf(__m256 x); /* _mm256_sin_ps(x) */
__m256 _ZGVcN8vv_powf(__m256 x, __m256 y); /* _mm256_pow_ps(x, y) */
void _ZGVcN8vvv_sincosf(__m256 x, __m256i ptrs_lo, __m256i ptrs_hi, __m256i ptrc_lo, __m256i ptrc_hi); /* _mm256_sincos_ps */
/* AVX2 */
__m256 _ZGVdN8v_cosf(__m256 x); /* _mm256_cos_ps(x) */
__m256 _ZGVdN8v_expf(__m256 x); /* _mm256_exp_ps(x) */
__m256 _ZGVdN8v_logf(__m256 x); /* _mm256_log_ps(x) */
__m256 _ZGVdN8v_sinf(__m256 x); /* _mm256_sin_ps(x) */
__m256 _ZGVdN8vv_powf(__m256 x, __m256 y); /* _mm256_pow_ps(x, y) */
void _ZGVdN8vvv_sincosf(__m256 x, __m256i ptrs_lo, __m256i ptrs_hi, __m256i ptrc_lo, __m256i ptrc_hi); /* _mm256_sincos_ps */
/* Note that I only tested the AVX2 function declarations. */
/* AVX-512: _ZGVeN16v_cosf, _ZGVeN8v_cos etc. */
有用的链接:
https://sourceware.org/glibc/wiki/libmvec
https://sourceware.org/glibc/wiki/libmvec?action=AttachFile&do=view&target=VectorABI.txt
https://github .com/sgallagher/glibc/blob/master/sysdeps/unix/sysv/linux/x86_64/libmvec.abilist
https://software.intel.com/sites/default/files/managed/b4/c8/Intel-Vector-Function-ABI.pdf
这篇关于SIMD寄存器的数学函数的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！
