为什么NaN值的比较与所有其他值的行为不同？也就是说,与运算符==,<=,> =,<,>的所有比较(其中一个或两个值为NaN)返回false,这与所有其他值的行为相反.
我想这可以通过某种方式简化数值计算,但我找不到明确说明的理由,甚至在Kahan 的IEEE 754状态讲义中也没有详细讨论其他设计决策.
这种异常行为在进行简单数据处理时会造成麻烦.例如,当在C程序中对某些实值字段的记录列表进行排序时,我需要编写额外的代码来处理NaN作为最大元素,否则排序算法可能会变得混乱.
编辑: 迄今为止的答案都认为比较NaNs毫无意义.
我同意,但这并不意味着正确的答案是错误的,而是一个非布尔值(NaB),幸运的是它不存在.
因此,在我看来,选择返回真或假的比较是任意的,对于一般数据处理,如果它遵循通常的定律(==的反射性,<= =,>的三分法),那将是有利的,以免数据结构依赖这些法律变得困惑.
因此,我要求打破这些法律的一些具体优势,而不仅仅是哲学推理.
编辑2: 我想我现在理解为什么使NaN最大化是一个坏主意,它会搞砸上限的计算.
可能需要NaN!= NaN以避免检测环路中的收敛,例如
while (x != oldX) {
oldX = x;
x = better_approximation(x);
}
然而,最好通过比较绝对差异与小限制来编写.所以恕我直言,这是打破NaN反身性的一个相对弱的论据.
我已经了解到一些Intel/AMD CPU可以同时进行多次复用并添加SSE/AVX:
每个周期的FLOPS用于沙桥和haswell SSE2/AVX/AVX2.
我想知道如何在代码中做到最好,我也想知道它是如何在CPU内部完成的.我的意思是超标量架构.假设我想做一个很长的总和,如下面的SSE:
//sum = a1*b1 + a2*b2 + a3*b3 +... where a is a scalar and b is a SIMD vector (e.g. from matrix multiplication)
sum = _mm_set1_ps(0.0f);
a1 = _mm_set1_ps(a[0]);
b1 = _mm_load_ps(&b[0]);
sum = _mm_add_ps(sum, _mm_mul_ps(a1, b1));
a2 = _mm_set1_ps(a[1]);
b2 = _mm_load_ps(&b[4]);
sum = _mm_add_ps(sum, _mm_mul_ps(a2, b2));
a3 = _mm_set1_ps(a[2]);
b3 = _mm_load_ps(&b[8]);
sum = _mm_add_ps(sum, _mm_mul_ps(a3, b3));
...
我的问题是如何将其转换为同时乘法并添加？数据可以依赖吗？我的意思是CPU可以_mm_add_ps(sum, _mm_mul_ps(a1, b1))同时执行还是在乘法中使用的寄存器和add必须是独立的？
最后,这如何适用于FMA(与Haswell)？是_mm_add_ps(sum, _mm_mul_ps(a1, b1))自动转换为单个FMA指令还是微操作？
C99添加了一个宏__STDC_IEC_559__,可用于测试编译器和标准库是否符合ISO/IEC/IEEE 60559(或IEEE 754)标准.
根据这个问题的答案,
如何检查-ieee-754-单精度-32位浮点表示,大多数C编译器不设置预处理器宏__STDC_IEC_559__.
根据海湾合作委员会的文件,它没有定义__STDC_IEC_559__.
我用GCC 4.9.2和Clang 3.6.0测试了这两个,使用glibc2.21使用以下代码.
//test.c
//#include <features.h>
int main(void) {
#if defined ( __STDC_IEC_559__ )
//#if defined ( __GCC_IEC_559__ )
return 1;
#else
return 0;
#endif
}
然后
echo $?
这表明此代码__STDC_IEC_559__是使用GCC定义的,而不是使用Clang定义的.然后我做了gcc -E,它显示该文件stdc-predef.h包含在内.这个文件定义__STDC_IEC_559__.
/* glibc's intent is to support the IEC 559 math functionality, real
and complex. If the GCC (4.9 and later) predefined macros
specifying compiler intent …我有一个i5-4250U,它有AVX2和FMA3.我正在测试Linux上的GCC 4.8.1中的一些密集矩阵乘法代码.下面是我编译的三种不同方式的列表.
SSE2: gcc matrix.cpp -o matrix_gcc -O3 -msse2 -fopenmp
AVX: gcc matrix.cpp -o matrix_gcc -O3 -mavx -fopenmp
AVX2+FMA: gcc matrix.cpp -o matrix_gcc -O3 -march=native -fopenmp -ffast-math
SSE2和AVX版本的性能明显不同.但是,AVX2 + FMA并不比AVX版本好.我不明白这一点.假设没有FMA,我获得了超过80%的CPU峰值触发器,但我认为我应该能够用FMA做得更好.矩阵乘法应直接受益于FMA.我基本上是在AVX中同时做八个点产品.当我检查march=native它给出:
cc -march=native -E -v - </dev/null 2>&1 | grep cc1 | grep fma
...-march=core-avx2 -mavx -mavx2 -mfma -mno-fma4 -msse4.2 -msse4.1 ...
所以我可以看到它已启用(只是为了确保我添加-mfma但它没有区别). ffast-math应该允许宽松的浮点模型如何在SSE/AVX中使用融合乘法 - 加法(FMA)指令
编辑:
基于Mysticial的评论我继续使用_mm256_fmadd_ps,现在AVX2 + FMA版本更快. 我不确定为什么编译器不会为我这样做. 对于超过1000x1000的矩阵,我现在得到大约80 GFLOPS(没有FMA的110%的峰值触发器).如果有人不信任我的峰值翻牌计算,这就是我所做的.
peak flops (no FMA) = frequency * simd_width * ILP * cores …c++当使用-O0或编译时，以下程序给出的数值不同-O1。
#include <iostream>
#include <array>
#include <cmath>
#include <iomanip>
int main()
{
std::array<double, 2> v;
v.fill(0);
std::get<0>(v) = 0x1.5b4d3afe1f7d2p-1;
double theta = 0x1.1aef12f82caf9p+2;
double c = std::cos(theta);
double s = std::sin(theta);
double vi = c * std::get<0>(v) - s * std::get<1>(v);
double vj = s * std::get<0>(v) + c * std::get<1>(v);
std::get<0>(v) = vi;
std::get<1>(v) = vj;
std::cout << std::hexfloat;
for (const auto& x : v)
std::cout << x << " ";
std::cout << std::endl; …根据文档,有一个fma()功能math.h.这非常好,我知道FMA如何工作以及如何使用它.但是,我不太确定这在实践中如何实施？我最感兴趣的是x86和x86_64架构.
是否存在FMA的浮点(非向量)指令,可能是IEEE-754 2008定义的？
是使用FMA3还是FMA4指令？
在依赖精度的情况下,是否存在确保使用真实FMA的内在因素？
在英特尔架构指令集扩展编程参考的2.5.3节"广播"中我们学习的比AVX512(和骑士角落)有
用于编码一些加载操作指令的数据广播的位字段,即从存储器加载数据并执行一些计算或数据移动操作的指令.
例如,使用英特尔汇编语法,我们可以在存储的地址广播标量,rax然后乘以16浮点数zmm2并将结果写入zmm1这样的
vmulps zmm1, zmm2, [rax] {1to16}
但是,没有内在的东西可以做到这一点.因此,使用内在函数,编译器应该能够折叠
__m512 bb = _mm512_set1_ps(b);
__m512 ab = _mm512_mul_ps(a,bb);
一个指令
vmulps zmm1, zmm2, [rax] {1to16}
但我没有观察到GCC这样做.我发现了一个GCC错误报告.
我观察到与GCC类似的FMA.例如,GCC 4.9不会崩溃_mm256_add_ps(_mm256_mul_ps(areg0,breg0) 为单个fma指令-Ofast.但是,GCC 5.1确实将它崩溃为一个单一的fma.至少有内在函数可以用FMA做到这一点,例如_mm256_fmadd_ps.但是没有例如_mm512_mulbroad_ps(vector,scalar)内在的.
海湾合作委员会可能会在某个时候解决这个问题,但在此之前,装配是唯
所以我的问题是如何在GCC中进行内联汇编？
我想我可能已经为上面的例子提出了GCC内联汇编的正确语法(但我不确定).
"vmulps (%%rax)%{1to16}, %%zmm1, %%zmm2\n\t"
我真的在寻找这样的功能
static inline __m512 mul_broad(__m512 a, float b) {
return a*b;
}
如果b在内存中指向rax它产生
vmulps (%rax){1to16}, %zmm0, %zmm0
ret
如果b它在xmm1它产生
vbroadcastss …我对-ffp-contractGNU GCC中的标志有疑问（请参阅https://gcc.gnu.org/onlinedocs/gcc/Optimize-Options.html）。
标志文档编写如下：
-ffp-contract=off禁用浮点表达式收缩。-ffp-contract=fast如果目标对它们有本机支持，则启用浮点表达式收缩，例如形成融合的乘法加法运算。-ffp-contract=on如果语言标准允许，则启用浮点表达式收缩。目前尚未实现，并且将其视为-ffp-contract=off。默认值为-ffp-contract=fast。
现在的问题是：
使用Haswell的FMA指令考虑以下指令序列:
__m256 r1 = _mm256_xor_ps (r1, r1);
r1 = _mm256_fmadd_ps (rp1, m6, r1);
r1 = _mm256_fmadd_ps (rp2, m7, r1);
r1 = _mm256_fmadd_ps (rp3, m8, r1);
__m256 r2 = _mm256_xor_ps (r2, r2);
r2 = _mm256_fmadd_ps (rp1, m3, r2);
r2 = _mm256_fmadd_ps (rp2, m4, r2);
r2 = _mm256_fmadd_ps (rp3, m5, r2);
__m256 r3 = _mm256_xor_ps (r3, r3);
r3 = _mm256_fmadd_ps (rp1, m0, r3);
r3 = _mm256_fmadd_ps (rp2, m1, r3);
r3 = _mm256_fmadd_ps (rp3, m2, r3);
可以使用非FMA指令表达相同的计算,如下所示:
__m256 i1 = _mm256_mul_ps (rp1, m6); …