有一种英特尔指令叫做F2XM1
,用于计算-1.0<x<1.0
时的2^x-1
,可用于计算tanh
,至少在某些范围内。然而,如果Agner's tables是可信的,这个操作的成本约为60 FLOPs。
另一个问题是向量化 - 据我所知,普通的glibc实现没有进行向量化。因此,如果您的程序使用了向量化,并且必须使用未向量化的tanh
实现,则会使程序变得更慢。为此,英特尔编译器具有mkl库,其中对其他函数一样,对tanh
进行了向量化。
从表中可以看出,最大成本约为每个操作10个时钟周期(一个浮点操作的成本约为1个时钟周期)。
我猜你可以通过使用-ffast-math
编译器选项赢得一些FLOPS,这会导致程序更快但不太精确(这是Cuda或c / c ++的选项,不确定是否可用于python / numpy)。
这段C++代码用于生成图表数据(使用g++ -std=c++11 -O2编译)。它的目的不是给出精确数字,而是让人们对成本有第一印象。
#include <chrono>
#include <iostream>
#include <vector>
#include <math.h>
int main(){
const std::vector<double> starts={1e-30, 1e-18, 1e-16, 1e-10, 1e-5, 1e-2, 1e-1, 0.5, 0.7, 0.9, 1.0, 2.0, 10, 20, 23, 100,1e3, 1e4};
const double FACTOR=1.0+1e-11;
const size_t ITER=100000000;
//warm-up:
double res=1.0;
for(size_t i=0;i<4*ITER;i++){
res*=FACTOR;
}
//overhead:
auto begin = std::chrono::high_resolution_clock::now();
for(size_t i=0;i<ITER;i++){
res*=FACTOR;
}
auto end = std::chrono::high_resolution_clock::now();
auto overhead=std::chrono::duration_cast<std::chrono::nanoseconds>(end-begin).count();
//std::cout<<"overhead: "<<overhead<<"\n";
//experiments:
for(auto start : starts){
begin=std::chrono::high_resolution_clock::now();
for(size_t i=0;i<ITER;i++){
res*=tanh(start);
start*=FACTOR;
}
auto end = std::chrono::high_resolution_clock::now();
auto time_needed=std::chrono::duration_cast<std::chrono::nanoseconds>(end-begin).count();
std::cout<<start<<" "<<time_needed/overhead<<"\n";
}
//overhead check:
begin = std::chrono::high_resolution_clock::now();
for(size_t i=0;i<ITER;i++){
res*=FACTOR;
}
end = std::chrono::high_resolution_clock::now();
auto overhead_new=std::chrono::duration_cast<std::chrono::nanoseconds>(end-begin).count();
std::cerr<<"overhead check: "<<overhead/overhead_new<<"\n";
std::cerr<<res;//don't optimize anything out...
}
Interval 0 x_small x_medium x_large
tanh(x) | x | polynomial approx. | 1-(2/(1+exp(2x))) | 1
现在,这是一份1993年的报告,但我认为这里并没有太多改变。
binary32
格式。提问者表示NVIDIA GPU是一个感兴趣的平台。我将专注于使用CUDA的这些GPU,因为我不熟悉CUDA的Python绑定。从图灵架构(计算能力7.5)开始,GPU包括指令MUFU.TANH
,用于计算具有约16位精度的单精度双曲正切。多功能单元(MUFU)支持的单精度函数通常通过在存储在ROM中的表中进行二次插值来计算。据我所知,MUFU.TANH
在虚拟汇编语言PTX的级别上公开,但(截至CUDA 11.2)不作为设备函数内部函数。
但是,鉴于该功能在PTX级别上公开,我们可以轻松地使用一行内联汇编代码创建自己的内部函数:
// Compute hyperbolic tangent for >= sm75. maxulperr = 133.95290, maxrelerr = 1.1126e-5
__forceinline__ __device__ float __tanhf (float a)
{
asm ("tanh.approx.f32 %0,%1; \n\t" : "=f"(a) : "f"(a));
return a;
}
// like copysignf(); when first argument is known to be positive
__forceinline__ __device__ float copysignf_pos (float a, float b)
{
return __int_as_float (__float_as_int (a) | (__float_as_int (b) & 0x80000000));
}
// Compute hyperbolic tangent for < sm_75. maxulperr = 108.82848, maxrelerr = 9.3450e-6
__forceinline__ __device__ float __tanhf (float a)
{
const float L2E = 1.442695041f;
float e, r, s, t, d;
s = fabsf (a);
t = -L2E * 2.0f * s;
asm ("ex2.approx.ftz.f32 %0,%1;\n\t" : "=f"(e) : "f"(t));
d = e + 1.0f;
asm ("rcp.approx.ftz.f32 %0,%1;\n\t" : "=f"(r) : "f"(d));
r = fmaf (e, -r, r);
if (s < 4.997253418e-3f) r = a;
if (!isnan (a)) r = copysignf_pos (r, a);
return r;
}
sm_70
,然后使用cuobjdump --dump-sass
反汇编二进制文件,显示八个浮点指令。我们还可以看到生成的机器代码(SASS)是无分支的。MUFU.EX2
和MUFU.RCP
。超过一定的参数范围,结果将为±1。// Compute hyperbolic tangent. maxulperr = 1.81484, maxrelerr = 1.9547e-7
__forceinline__ __device__ float my_tanhf (float a)
{
const float L2E = 1.442695041f;
float p, s, t, r;
t = fabsf (a);
if (t >= 307.0f/512.0f) { // 0.599609375
r = L2E * 2.0f * t;
asm ("ex2.approx.ftz.f32 %0,%1;\n\t" : "=f"(r) : "f"(r));
r = 1.0f + r;
asm ("rcp.approx.ftz.f32 %0,%1;\n\t" : "=f"(r) : "f"(r));
r = fmaf (r, -2.0f, 1.0f);
if (t >= 9.03125f) r = 1.0f;
r = copysignf_pos (r, a);
} else {
s = a * a;
p = 1.57394409e-2f; // 0x1.01e000p-6
p = fmaf (p, s, -5.23025580e-2f); // -0x1.ac766ap-5
p = fmaf (p, s, 1.33152470e-1f); // 0x1.10b23ep-3
p = fmaf (p, s, -3.33327681e-1f); // -0x1.5553dap-2
p = fmaf (p, s, 0.0f);
r = fmaf (p, a, a);
}
return r;
}
sm75
目标机器码可以看出该分支被保留。这意味着在一般情况下,所有活动线程中的一些将遵循分支的一侧,而其余线程将遵循分支的另一侧,需要后续同步。因此,为了获得所需的计算工作量的实际概念,我们需要结合两个执行路径的浮点指令计数。这总共有十三条浮点指令。