是否有一个 C++ 函数 returns 正好是内置 CPU 运算 RSQRTSS 的平方根倒数的值？
Is there a C++ function that returns exactly the value of the built-in CPU operation RSQRTSS for inverse square root?
我正在寻找一个 C++ 函数，它 returns 一个浮点数的逆平方根：rsqrt(x) = 1/sqrt(x) 通过使用像内置 XMM 操作 RSQRTSS 这样的精确方法（cf. https://www.felixcloutier.com/x86/rsqrtss)。 （即，我想要内置的近似值而不是更精确的 1/sqrtf，而且我不关心速度（很多）。）
根据这个问题：
Is there a fast C or C++ standard library function for double precision inverse square root?
...至少没有 "fast way with double precision" 可以使用标准 C++ 库来完成此操作。但是如何做到慢，非标准和 floats？
标准库中没有执行此操作的函数，但您的编译器可能会优化表达式 1 / sqrt(value)，使其发出 RSQRTSS 指令。例如，使用编译器标志 -ffast-math -march=native，GCC 将发出该指令，请参阅：https://godbolt.org/z/cL6seG
就其价值而言，正如@François Andrieux 所建议的那样，我最终在 C++ 中以普通汇编实现了它（更准确地说，我使用了 ASMJIT）。
这很好用，尽管它有失去可移植性的缺点（虽然不如普通 asm）。但这在某种程度上是我的问题所固有的，因为我想使用一个非常具体的 x86 函数。
这是我的代码：
typedef float(*JITFunc)();
JITFunc func;
asmjit::JitRuntime jit_runtime;
asmjit::CodeHolder code;
code.init(jit_runtime.getCodeInfo());
asmjit::X86Compiler cc(&code);
cc.addFunc(asmjit::FuncSignature0<float>());
float value = 2.71; // Some example value.
asmjit::X86Xmm x = cc.newXmm();
setXmmVar(cc, x, value);
cc.rsqrtss(x, x); // THE asm function.
cc.ret(x);
cc.endFunc();
cc.finalize();
jit_runtime.add(&func, &code);
return func(); // Or something to that effect. func() is the result, anyway.
RSQRTSS 指令可通过 immintrin.h 中声明的 _mm_rsqrt_ss() 内在函数轻松访问。但是我们也可以在软件中模拟指令，就像下面的 my_rsqrtf() 函数所做的那样。通过简单地观察 RSQRTSS 的输出，很容易发现它的函数值基于 211 个条目的（虚拟）table，每个条目的大小为 12 位.
注意 "virtual" 属性，因为硬件不太可能使用直接 24 Kbit table。我对 table 条目中模式的分析不建议使用 bipartite tables。一种更简单的压缩方案——就像我在下面的代码中使用的那样——基于 table 基值和 table 偏移量 可能 是用过的。我的方案只需要一个窄加器，但将 ROM 存储减少到 13 Kbit，即几乎一半。
下面的实现是在使用 Ivy Bridge 架构的英特尔至强处理器 E3-1270 V2 上开发和测试的。 RSQRTSS 在各种 Intel 架构之间的实现可能存在一些功能差异，并且这种差异可能存在于来自不同 x86-84 供应商的架构之间。
下面的框架检查 my_rsqrtf() 的仿真为所有四种舍入模式、两种 DAZ（非正规数为零）模式和两种 FTZ（齐平）提供与 RSQRTSS 相同的结果到零）模式。我们发现函数结果不受任何模式的影响，这与英特尔在 RSQRTSS 中指定的方式相符
英特尔® 64 位和 IA-32 架构软件开发人员手册：
The RSQRTSS instruction is not affected by the rounding control bits in the MXCSR register. When a source value is a 0.0, an ∞ of the sign of the source value is returned. A denormal source value is treated as a 0.0 (of the same sign). When a source value is a negative value (other than −0.0), a floating-point indefinite is returned. When a source value is an SNaN or QNaN, the SNaN is converted to a QNaN or the source QNaN is returned.
#include <stdio.h>
#include <stdlib.h>
#include <stdint.h>
#include <string.h>
#include <math.h>
#include "immintrin.h"
/* SSE reference for RSQRTSS instruction */
float sse_rsqrtf (float a, uint32_t daz, uint32_t ftz, uint32_t rnd)
{
__m128 b, t;
float res;
uint32_t old_mxcsr;
old_mxcsr = _mm_getcsr();
_MM_SET_DENORMALS_ZERO_MODE (daz);
_MM_SET_FLUSH_ZERO_MODE (ftz);
_MM_SET_ROUNDING_MODE (rnd);
b = _mm_set_ss (a);
t = _mm_rsqrt_ss (b);
_mm_store_ss (&res, t);
_mm_setcsr (old_mxcsr);
return res;
}
inline uint32_t float_as_uint32 (float a)
{
uint32_t r;
memcpy (&r, &a, sizeof r);
return r;
}
inline float uint32_as_float (uint32_t a)
{
float r;
memcpy (&r, &a, sizeof r);
return r;
}
#define LOG2_NBR_TAB_ENTRIES (11)
#define NBR_TAB_ENTRIES (1 << LOG2_NBR_TAB_ENTRIES)
#define TAB_ENTRY_BITS (12)
#define BASE_TAB_ENTRY_BITS (9)
/* 128 9-bit entries = 1152 bits */
const uint16_t base_tab[128] = {
0x0ce, 0x0c9, 0x0c3, 0x0be, 0x0b9, 0x0b4, 0x0af, 0x0aa,
0x0a6, 0x0a1, 0x09d, 0x098, 0x094, 0x090, 0x08b, 0x087,
0x083, 0x07f, 0x07c, 0x078, 0x074, 0x070, 0x06d, 0x069,
0x066, 0x062, 0x05f, 0x05c, 0x058, 0x055, 0x052, 0x04f,
0x04c, 0x049, 0x046, 0x043, 0x040, 0x03d, 0x03a, 0x038,
0x035, 0x032, 0x030, 0x02d, 0x02a, 0x028, 0x025, 0x023,
0x021, 0x01e, 0x01c, 0x019, 0x017, 0x015, 0x013, 0x010,
0x00e, 0x00c, 0x00a, 0x008, 0x006, 0x004, 0x002, 0x000,
0x1f8, 0x1f0, 0x1e9, 0x1e1, 0x1da, 0x1d3, 0x1cc, 0x1c5,
0x1bf, 0x1b8, 0x1b2, 0x1ab, 0x1a5, 0x19f, 0x199, 0x194,
0x18e, 0x188, 0x183, 0x17e, 0x178, 0x173, 0x16e, 0x169,
0x164, 0x15f, 0x15a, 0x156, 0x151, 0x14d, 0x148, 0x144,
0x13f, 0x13b, 0x137, 0x133, 0x12f, 0x12b, 0x127, 0x123,
0x11f, 0x11b, 0x118, 0x114, 0x110, 0x10d, 0x109, 0x106,
0x102, 0x0ff, 0x0fc, 0x0f8, 0x0f5, 0x0f2, 0x0ef, 0x0eb,
0x0e8, 0x0e5, 0x0e2, 0x0df, 0x0dc, 0x0d9, 0x0d7, 0x0d4,
};
/* 2048 6-bit entries = 12288 bits */
const uint8_t ofs_tab[2048] = {
0x2f, 0x2c, 0x2a, 0x27, 0x24, 0x21, 0x1e, 0x1c, 0x19, 0x16, 0x13, 0x10, 0x0e, 0x0b, 0x08, 0x05,
0x2b, 0x28, 0x25, 0x22, 0x1f, 0x1d, 0x1a, 0x17, 0x15, 0x12, 0x0f, 0x0c, 0x0a, 0x07, 0x04, 0x02,
0x2f, 0x2c, 0x29, 0x27, 0x24, 0x21, 0x1f, 0x1c, 0x19, 0x17, 0x14, 0x11, 0x0f, 0x0c, 0x09, 0x07,
0x2c, 0x29, 0x27, 0x24, 0x22, 0x1f, 0x1c, 0x1a, 0x17, 0x15, 0x12, 0x0f, 0x0d, 0x0a, 0x08, 0x05,
0x2a, 0x28, 0x25, 0x23, 0x20, 0x1e, 0x1b, 0x18, 0x16, 0x13, 0x11, 0x0e, 0x0c, 0x09, 0x07, 0x04,
0x2a, 0x27, 0x24, 0x22, 0x1f, 0x1d, 0x1a, 0x18, 0x15, 0x13, 0x10, 0x0e, 0x0b, 0x09, 0x07, 0x04,
0x2a, 0x27, 0x25, 0x22, 0x20, 0x1d, 0x1b, 0x18, 0x16, 0x13, 0x11, 0x0f, 0x0c, 0x0a, 0x07, 0x05,
0x2a, 0x28, 0x26, 0x23, 0x21, 0x1e, 0x1c, 0x1a, 0x17, 0x15, 0x12, 0x10, 0x0e, 0x0b, 0x09, 0x07,
0x24, 0x22, 0x1f, 0x1d, 0x1b, 0x18, 0x16, 0x14, 0x11, 0x0f, 0x0d, 0x0a, 0x08, 0x06, 0x03, 0x01,
0x27, 0x24, 0x22, 0x20, 0x1d, 0x1b, 0x19, 0x16, 0x14, 0x12, 0x10, 0x0d, 0x0b, 0x09, 0x06, 0x04,
0x22, 0x20, 0x1d, 0x1b, 0x19, 0x17, 0x14, 0x12, 0x10, 0x0e, 0x0b, 0x09, 0x07, 0x05, 0x02, 0x00,
0x26, 0x24, 0x21, 0x1f, 0x1d, 0x1b, 0x19, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0b, 0x09, 0x07, 0x05,
0x23, 0x20, 0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x13, 0x11, 0x0f, 0x0d, 0x0b, 0x09, 0x06, 0x04, 0x02,
0x20, 0x1e, 0x1c, 0x1a, 0x17, 0x15, 0x13, 0x11, 0x0f, 0x0d, 0x0b, 0x09, 0x06, 0x04, 0x02, 0x00,
0x26, 0x24, 0x22, 0x20, 0x1e, 0x1c, 0x19, 0x17, 0x15, 0x13, 0x11, 0x0f, 0x0d, 0x0b, 0x09, 0x07,
0x25, 0x23, 0x21, 0x1f, 0x1d, 0x1a, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x08, 0x06,
0x24, 0x22, 0x20, 0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x08, 0x06,
0x24, 0x22, 0x20, 0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x08, 0x06,
0x1d, 0x1b, 0x19, 0x17, 0x15, 0x13, 0x11, 0x0f, 0x0d, 0x0b, 0x09, 0x07, 0x05, 0x03, 0x01, 0x00,
0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0b, 0x09, 0x07, 0x05, 0x03, 0x01,
0x1f, 0x1d, 0x1c, 0x1a, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0d, 0x0b, 0x09, 0x07, 0x05, 0x03,
0x21, 0x20, 0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x15, 0x13, 0x11, 0x0f, 0x0d, 0x0b, 0x0a, 0x08, 0x06,
0x1c, 0x1a, 0x19, 0x17, 0x15, 0x13, 0x11, 0x10, 0x0e, 0x0c, 0x0a, 0x08, 0x07, 0x05, 0x03, 0x01,
0x1f, 0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x15, 0x13, 0x11, 0x0f, 0x0e, 0x0c, 0x0a, 0x08, 0x07, 0x05,
0x1b, 0x19, 0x18, 0x16, 0x14, 0x12, 0x11, 0x0f, 0x0d, 0x0b, 0x0a, 0x08, 0x06, 0x04, 0x03, 0x01,
0x1f, 0x1e, 0x1c, 0x1a, 0x18, 0x17, 0x15, 0x13, 0x12, 0x10, 0x0e, 0x0c, 0x0b, 0x09, 0x07, 0x06,
0x1c, 0x1a, 0x19, 0x17, 0x15, 0x13, 0x12, 0x10, 0x0e, 0x0d, 0x0b, 0x09, 0x08, 0x06, 0x04, 0x03,
0x19, 0x17, 0x16, 0x14, 0x12, 0x11, 0x0f, 0x0d, 0x0c, 0x0a, 0x08, 0x07, 0x05, 0x03, 0x02, 0x00,
0x1f, 0x1d, 0x1b, 0x1a, 0x18, 0x16, 0x15, 0x13, 0x11, 0x10, 0x0e, 0x0d, 0x0b, 0x09, 0x08, 0x06,
0x1d, 0x1b, 0x19, 0x18, 0x16, 0x14, 0x13, 0x11, 0x10, 0x0e, 0x0c, 0x0b, 0x09, 0x08, 0x06, 0x04,
0x1b, 0x19, 0x18, 0x16, 0x15, 0x13, 0x11, 0x10, 0x0e, 0x0d, 0x0b, 0x0a, 0x08, 0x06, 0x05, 0x03,
0x1a, 0x18, 0x17, 0x15, 0x13, 0x12, 0x10, 0x0f, 0x0d, 0x0c, 0x0a, 0x09, 0x07, 0x06, 0x04, 0x02,
0x19, 0x17, 0x16, 0x14, 0x13, 0x11, 0x10, 0x0e, 0x0d, 0x0b, 0x0a, 0x08, 0x07, 0x05, 0x03, 0x02,
0x18, 0x17, 0x15, 0x14, 0x12, 0x11, 0x0f, 0x0e, 0x0c, 0x0b, 0x09, 0x08, 0x06, 0x05, 0x03, 0x02,
0x18, 0x17, 0x15, 0x14, 0x12, 0x11, 0x0f, 0x0e, 0x0d, 0x0b, 0x0a, 0x08, 0x07, 0x05, 0x04, 0x02,
0x19, 0x17, 0x16, 0x14, 0x13, 0x11, 0x10, 0x0e, 0x0d, 0x0c, 0x0a, 0x09, 0x07, 0x06, 0x04, 0x03,
0x19, 0x18, 0x16, 0x15, 0x14, 0x12, 0x11, 0x0f, 0x0e, 0x0c, 0x0b, 0x0a, 0x08, 0x07, 0x05, 0x04,
0x1a, 0x19, 0x18, 0x16, 0x15, 0x13, 0x12, 0x10, 0x0f, 0x0e, 0x0c, 0x0b, 0x09, 0x08, 0x07, 0x05,
0x1c, 0x1a, 0x19, 0x18, 0x16, 0x15, 0x13, 0x12, 0x11, 0x0f, 0x0e, 0x0c, 0x0b, 0x0a, 0x08, 0x07,
0x15, 0x14, 0x13, 0x11, 0x10, 0x0f, 0x0d, 0x0c, 0x0a, 0x09, 0x08, 0x06, 0x05, 0x04, 0x02, 0x01,
0x17, 0x16, 0x15, 0x13, 0x12, 0x11, 0x0f, 0x0e, 0x0d, 0x0b, 0x0a, 0x08, 0x07, 0x06, 0x04, 0x03,
0x1a, 0x18, 0x17, 0x16, 0x14, 0x13, 0x12, 0x10, 0x0f, 0x0e, 0x0c, 0x0b, 0x0a, 0x08, 0x07, 0x06,
0x14, 0x13, 0x12, 0x10, 0x0f, 0x0e, 0x0c, 0x0b, 0x0a, 0x08, 0x07, 0x06, 0x05, 0x03, 0x02, 0x01,
0x17, 0x16, 0x15, 0x13, 0x12, 0x11, 0x0f, 0x0e, 0x0d, 0x0c, 0x0a, 0x09, 0x08, 0x06, 0x05, 0x04,
0x1b, 0x19, 0x18, 0x17, 0x15, 0x14, 0x13, 0x12, 0x10, 0x0f, 0x0e, 0x0c, 0x0b, 0x0a, 0x09, 0x07,
0x16, 0x15, 0x13, 0x12, 0x11, 0x10, 0x0e, 0x0d, 0x0c, 0x0b, 0x09, 0x08, 0x07, 0x06, 0x04, 0x03,
0x1a, 0x19, 0x17, 0x16, 0x15, 0x14, 0x12, 0x11, 0x10, 0x0f, 0x0d, 0x0c, 0x0b, 0x0a, 0x08, 0x07,
0x16, 0x15, 0x13, 0x12, 0x11, 0x10, 0x0e, 0x0d, 0x0c, 0x0b, 0x0a, 0x08, 0x07, 0x06, 0x05, 0x03,
0x12, 0x11, 0x10, 0x0f, 0x0d, 0x0c, 0x0b, 0x0a, 0x08, 0x07, 0x06, 0x05, 0x04, 0x02, 0x01, 0x00,
0x17, 0x16, 0x14, 0x13, 0x12, 0x11, 0x10, 0x0e, 0x0d, 0x0c, 0x0b, 0x0a, 0x08, 0x07, 0x06, 0x05,
0x14, 0x12, 0x11, 0x10, 0x0f, 0x0e, 0x0d, 0x0b, 0x0a, 0x09, 0x08, 0x07, 0x05, 0x04, 0x03, 0x02,
0x19, 0x18, 0x16, 0x15, 0x14, 0x13, 0x12, 0x11, 0x0f, 0x0e, 0x0d, 0x0c, 0x0b, 0x0a, 0x08, 0x07,
0x16, 0x15, 0x14, 0x13, 0x11, 0x10, 0x0f, 0x0e, 0x0d, 0x0c, 0x0b, 0x09, 0x08, 0x07, 0x06, 0x05,
0x14, 0x13, 0x11, 0x10, 0x0f, 0x0e, 0x0d, 0x0c, 0x0b, 0x09, 0x08, 0x07, 0x06, 0x05, 0x04, 0x03,
0x11, 0x10, 0x0f, 0x0e, 0x0d, 0x0c, 0x0b, 0x0a, 0x08, 0x07, 0x06, 0x05, 0x04, 0x03, 0x02, 0x01,
0x18, 0x16, 0x15, 0x14, 0x13, 0x12, 0x11, 0x10, 0x0f, 0x0e, 0x0c, 0x0b, 0x0a, 0x09, 0x08, 0x07,
0x16, 0x15, 0x14, 0x12, 0x11, 0x10, 0x0f, 0x0e, 0x0d, 0x0c, 0x0b, 0x0a, 0x09, 0x08, 0x06, 0x05,
0x14, 0x13, 0x12, 0x11, 0x10, 0x0f, 0x0e, 0x0d, 0x0c, 0x0b, 0x09, 0x08, 0x07, 0x06, 0x05, 0x04,
0x13, 0x12, 0x11, 0x10, 0x0f, 0x0e, 0x0d, 0x0b, 0x0a, 0x09, 0x08, 0x07, 0x06, 0x05, 0x04, 0x03,
0x12, 0x11, 0x10, 0x0f, 0x0e, 0x0d, 0x0c, 0x0a, 0x09, 0x08, 0x07, 0x06, 0x05, 0x04, 0x03, 0x02,
0x11, 0x10, 0x0f, 0x0e, 0x0d, 0x0c, 0x0b, 0x0a, 0x09, 0x08, 0x07, 0x06, 0x04, 0x03, 0x02, 0x01,
0x10, 0x0f, 0x0e, 0x0d, 0x0c, 0x0b, 0x0a, 0x09, 0x08, 0x07, 0x06, 0x05, 0x04, 0x03, 0x02, 0x01,
0x10, 0x0f, 0x0e, 0x0d, 0x0c, 0x0b, 0x0a, 0x09, 0x08, 0x07, 0x06, 0x05, 0x04, 0x03, 0x02, 0x01,
0x10, 0x0f, 0x0e, 0x0d, 0x0c, 0x0b, 0x0a, 0x09, 0x08, 0x07, 0x06, 0x05, 0x04, 0x03, 0x02, 0x01,
0x3e, 0x3a, 0x36, 0x32, 0x2e, 0x2a, 0x26, 0x22, 0x1e, 0x1a, 0x16, 0x12, 0x0e, 0x0b, 0x07, 0x03,
0x3f, 0x3b, 0x37, 0x33, 0x2f, 0x2b, 0x27, 0x24, 0x20, 0x1c, 0x18, 0x14, 0x10, 0x0c, 0x09, 0x05,
0x39, 0x35, 0x31, 0x2e, 0x2a, 0x26, 0x22, 0x1e, 0x1b, 0x17, 0x13, 0x0f, 0x0c, 0x08, 0x04, 0x00,
0x3d, 0x39, 0x35, 0x31, 0x2e, 0x2a, 0x26, 0x23, 0x1f, 0x1b, 0x18, 0x14, 0x10, 0x0d, 0x09, 0x05,
0x3a, 0x36, 0x32, 0x2f, 0x2b, 0x27, 0x24, 0x20, 0x1d, 0x19, 0x15, 0x12, 0x0e, 0x0b, 0x07, 0x03,
0x38, 0x34, 0x31, 0x2d, 0x2a, 0x26, 0x22, 0x1f, 0x1b, 0x18, 0x14, 0x11, 0x0d, 0x0a, 0x06, 0x03,
0x37, 0x34, 0x30, 0x2d, 0x29, 0x26, 0x22, 0x1f, 0x1b, 0x18, 0x15, 0x11, 0x0e, 0x0a, 0x07, 0x03,
0x38, 0x35, 0x31, 0x2e, 0x2a, 0x27, 0x24, 0x20, 0x1d, 0x19, 0x16, 0x13, 0x0f, 0x0c, 0x09, 0x05,
0x32, 0x2e, 0x2b, 0x28, 0x24, 0x21, 0x1e, 0x1a, 0x17, 0x14, 0x11, 0x0d, 0x0a, 0x07, 0x03, 0x00,
0x35, 0x31, 0x2e, 0x2b, 0x28, 0x24, 0x21, 0x1e, 0x1b, 0x17, 0x14, 0x11, 0x0e, 0x0a, 0x07, 0x04,
0x31, 0x2e, 0x2a, 0x27, 0x24, 0x21, 0x1e, 0x1a, 0x17, 0x14, 0x11, 0x0e, 0x0b, 0x07, 0x04, 0x01,
0x36, 0x33, 0x30, 0x2c, 0x29, 0x26, 0x23, 0x20, 0x1d, 0x1a, 0x17, 0x13, 0x10, 0x0d, 0x0a, 0x07,
0x34, 0x31, 0x2e, 0x2b, 0x28, 0x25, 0x21, 0x1e, 0x1b, 0x18, 0x15, 0x12, 0x0f, 0x0c, 0x09, 0x06,
0x33, 0x30, 0x2d, 0x2a, 0x27, 0x24, 0x21, 0x1e, 0x1b, 0x18, 0x15, 0x12, 0x0f, 0x0c, 0x09, 0x06,
0x33, 0x30, 0x2d, 0x2a, 0x27, 0x24, 0x21, 0x1e, 0x1b, 0x18, 0x15, 0x13, 0x10, 0x0d, 0x0a, 0x07,
0x2c, 0x29, 0x26, 0x23, 0x20, 0x1d, 0x1a, 0x18, 0x15, 0x12, 0x0f, 0x0c, 0x09, 0x06, 0x03, 0x01,
0x2e, 0x2b, 0x28, 0x25, 0x22, 0x1f, 0x1d, 0x1a, 0x17, 0x14, 0x11, 0x0e, 0x0c, 0x09, 0x06, 0x03,
0x30, 0x2e, 0x2b, 0x28, 0x25, 0x22, 0x20, 0x1d, 0x1a, 0x17, 0x14, 0x12, 0x0f, 0x0c, 0x09, 0x07,
0x2c, 0x29, 0x26, 0x24, 0x21, 0x1e, 0x1b, 0x19, 0x16, 0x13, 0x10, 0x0e, 0x0b, 0x08, 0x06, 0x03,
0x28, 0x25, 0x23, 0x20, 0x1d, 0x1b, 0x18, 0x15, 0x13, 0x10, 0x0d, 0x0b, 0x08, 0x05, 0x03, 0x00,
0x2d, 0x2b, 0x28, 0x25, 0x23, 0x20, 0x1d, 0x1b, 0x18, 0x15, 0x13, 0x10, 0x0e, 0x0b, 0x08, 0x06,
0x2b, 0x28, 0x26, 0x23, 0x21, 0x1e, 0x1b, 0x19, 0x16, 0x14, 0x11, 0x0f, 0x0c, 0x09, 0x07, 0x04,
0x2a, 0x27, 0x25, 0x22, 0x1f, 0x1d, 0x1a, 0x18, 0x15, 0x13, 0x10, 0x0e, 0x0b, 0x09, 0x06, 0x03,
0x29, 0x26, 0x24, 0x21, 0x1f, 0x1c, 0x1a, 0x17, 0x15, 0x12, 0x10, 0x0d, 0x0b, 0x08, 0x06, 0x03,
0x29, 0x26, 0x24, 0x21, 0x1f, 0x1d, 0x1a, 0x18, 0x15, 0x13, 0x10, 0x0e, 0x0b, 0x09, 0x06, 0x04,
0x2a, 0x27, 0x25, 0x22, 0x20, 0x1d, 0x1b, 0x19, 0x16, 0x14, 0x11, 0x0f, 0x0d, 0x0a, 0x08, 0x05,
0x2b, 0x29, 0x26, 0x24, 0x21, 0x1f, 0x1d, 0x1a, 0x18, 0x15, 0x13, 0x11, 0x0e, 0x0c, 0x0a, 0x07,
0x25, 0x23, 0x20, 0x1e, 0x1b, 0x19, 0x17, 0x14, 0x12, 0x10, 0x0d, 0x0b, 0x09, 0x06, 0x04, 0x02,
0x27, 0x25, 0x23, 0x20, 0x1e, 0x1c, 0x1a, 0x17, 0x15, 0x13, 0x10, 0x0e, 0x0c, 0x09, 0x07, 0x05,
0x23, 0x20, 0x1e, 0x1c, 0x1a, 0x17, 0x15, 0x13, 0x10, 0x0e, 0x0c, 0x0a, 0x07, 0x05, 0x03, 0x01,
0x26, 0x24, 0x22, 0x20, 0x1d, 0x1b, 0x19, 0x17, 0x15, 0x12, 0x10, 0x0e, 0x0c, 0x09, 0x07, 0x05,
0x23, 0x21, 0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x13, 0x11, 0x0f, 0x0d, 0x0b, 0x08, 0x06, 0x04, 0x02,
0x28, 0x25, 0x23, 0x21, 0x1f, 0x1d, 0x1b, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x09, 0x07,
0x25, 0x23, 0x21, 0x1f, 0x1d, 0x1a, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x09, 0x07, 0x05,
0x23, 0x21, 0x1f, 0x1d, 0x1b, 0x19, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x08, 0x06, 0x04,
0x22, 0x1f, 0x1d, 0x1b, 0x19, 0x17, 0x15, 0x13, 0x11, 0x0f, 0x0d, 0x0b, 0x09, 0x07, 0x05, 0x03,
0x21, 0x1f, 0x1c, 0x1a, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x08, 0x06, 0x04, 0x02,
0x20, 0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x08, 0x06, 0x04, 0x02,
0x20, 0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x08, 0x06, 0x04, 0x02,
0x20, 0x1f, 0x1d, 0x1b, 0x19, 0x17, 0x15, 0x13, 0x11, 0x0f, 0x0d, 0x0b, 0x09, 0x07, 0x05, 0x03,
0x21, 0x1f, 0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x08, 0x07, 0x05,
0x23, 0x21, 0x1f, 0x1d, 0x1b, 0x19, 0x17, 0x15, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x08, 0x06,
0x1c, 0x1b, 0x19, 0x17, 0x15, 0x13, 0x11, 0x0f, 0x0e, 0x0c, 0x0a, 0x08, 0x06, 0x04, 0x02, 0x01,
0x1f, 0x1d, 0x1b, 0x19, 0x17, 0x15, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x09, 0x07, 0x05, 0x03,
0x21, 0x1f, 0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x15, 0x13, 0x11, 0x0f, 0x0d, 0x0c, 0x0a, 0x08, 0x06,
0x1c, 0x1b, 0x19, 0x17, 0x15, 0x13, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x09, 0x07, 0x05, 0x03, 0x02,
0x20, 0x1e, 0x1c, 0x1a, 0x19, 0x17, 0x15, 0x13, 0x12, 0x10, 0x0e, 0x0c, 0x0b, 0x09, 0x07, 0x05,
0x1c, 0x1a, 0x18, 0x16, 0x15, 0x13, 0x11, 0x0f, 0x0e, 0x0c, 0x0a, 0x08, 0x07, 0x05, 0x03, 0x01,
0x20, 0x1e, 0x1c, 0x1b, 0x19, 0x17, 0x15, 0x14, 0x12, 0x10, 0x0f, 0x0d, 0x0b, 0x09, 0x08, 0x06,
0x1c, 0x1b, 0x19, 0x17, 0x15, 0x14, 0x12, 0x10, 0x0f, 0x0d, 0x0b, 0x0a, 0x08, 0x06, 0x05, 0x03,
0x19, 0x17, 0x16, 0x14, 0x12, 0x11, 0x0f, 0x0d, 0x0c, 0x0a, 0x08, 0x07, 0x05, 0x03, 0x02, 0x00,
0x1e, 0x1d, 0x1b, 0x19, 0x18, 0x16, 0x14, 0x13, 0x11, 0x10, 0x0e, 0x0c, 0x0b, 0x09, 0x07, 0x06,
0x1c, 0x1a, 0x19, 0x17, 0x16, 0x14, 0x12, 0x11, 0x0f, 0x0d, 0x0c, 0x0a, 0x08, 0x07, 0x05, 0x04,
0x1a, 0x18, 0x17, 0x15, 0x14, 0x12, 0x10, 0x0f, 0x0d, 0x0b, 0x0a, 0x08, 0x07, 0x05, 0x03, 0x02,
0x18, 0x17, 0x15, 0x13, 0x12, 0x10, 0x0f, 0x0d, 0x0c, 0x0a, 0x08, 0x07, 0x05, 0x04, 0x02, 0x00,
0x1f, 0x1d, 0x1c, 0x1a, 0x19, 0x17, 0x15, 0x14, 0x12, 0x11, 0x0f, 0x0e, 0x0c, 0x0a, 0x09, 0x07,
0x1e, 0x1c, 0x1b, 0x19, 0x18, 0x16, 0x14, 0x13, 0x11, 0x10, 0x0e, 0x0d, 0x0b, 0x0a, 0x08, 0x07,
0x1d, 0x1c, 0x1a, 0x18, 0x17, 0x15, 0x14, 0x12, 0x11, 0x0f, 0x0e, 0x0c, 0x0b, 0x09, 0x08, 0x06,
0x1d, 0x1b, 0x1a, 0x18, 0x17, 0x15, 0x13, 0x12, 0x10, 0x0f, 0x0d, 0x0c, 0x0a, 0x09, 0x07, 0x06,
0x1c, 0x1b, 0x19, 0x18, 0x16, 0x15, 0x13, 0x12, 0x10, 0x0f, 0x0d, 0x0c, 0x0b, 0x09, 0x08, 0x06,
0x1d, 0x1b, 0x1a, 0x18, 0x17, 0x15, 0x14, 0x12, 0x11, 0x0f, 0x0e, 0x0c, 0x0b, 0x09, 0x08, 0x06,
0x1d, 0x1c, 0x1a, 0x19, 0x17, 0x16, 0x14, 0x13, 0x11, 0x10, 0x0e, 0x0d, 0x0b, 0x0a, 0x09, 0x07,
0x16, 0x14, 0x13, 0x11, 0x10, 0x0e, 0x0d, 0x0c, 0x0a, 0x09, 0x07, 0x06, 0x04, 0x03, 0x02, 0x00,
0x17, 0x15, 0x14, 0x12, 0x11, 0x10, 0x0e, 0x0d, 0x0b, 0x0a, 0x08, 0x07, 0x06, 0x04, 0x03, 0x01,
};
#define IEEE_BINARY32_EXPO_BIAS (127)
#define IEEE_BINARY32_MANT_BITS (23)
#define IEEE_BINARY32_EXPO_BITS (8)
#define IEEE_BINARY32_EXPO_MASK (0x7f800000)
#define IEEE_BINARY32_NAN_INDEF (0xffc00000)
#define IEEE_BINARY32_POS_INF (0x7f800000)
#define IEEE_BINARY32_POS_ZERO (0x00000000)
#define IEEE_BINARY32_MIN_NORM (0x00800000)
#define IEEE_BINARY32_SIGN_BIT (0x80000000)
/* Emulate the RSQRTSS instruction in software */
float my_rsqrtf (float x)
{
float r;
uint32_t arg, res, idx, expo, mant;
arg = float_as_uint32 (x);
/* zeros and subnormals */
if ((arg & ~IEEE_BINARY32_SIGN_BIT) < IEEE_BINARY32_MIN_NORM) {
res = IEEE_BINARY32_POS_INF | (arg & IEEE_BINARY32_SIGN_BIT);
r = uint32_as_float (res);
}
/* NaNs */
else if ((arg & ~IEEE_BINARY32_SIGN_BIT) > IEEE_BINARY32_POS_INF) {
r = x + x; // convert SNaN to QNaN
}
/* negative arguments */
else if (arg & IEEE_BINARY32_SIGN_BIT) {
res = IEEE_BINARY32_NAN_INDEF;
r = uint32_as_float (res);
}
/* positive infinity */
else if (arg == IEEE_BINARY32_POS_INF) {
res = IEEE_BINARY32_POS_ZERO;
r = uint32_as_float (res);
}
/* positive normals */
else {
/* extract exponent lsb and leading mantissa bits for table index */
expo = (arg & IEEE_BINARY32_EXPO_MASK) >> IEEE_BINARY32_MANT_BITS;
idx = (arg >> (IEEE_BINARY32_MANT_BITS - LOG2_NBR_TAB_ENTRIES + 1))
& (NBR_TAB_ENTRIES - 1);
/* compute exponent and mantissa of reciprocal square root */
expo = (3 * IEEE_BINARY32_EXPO_BIAS + ~expo) >> 1;
mant = (((base_tab [idx >> 4] << (TAB_ENTRY_BITS - BASE_TAB_ENTRY_BITS)) + ofs_tab [idx])
<< (IEEE_BINARY32_MANT_BITS - TAB_ENTRY_BITS));
/* combine exponent and mantissa bits to compute final result */
res = (expo << IEEE_BINARY32_MANT_BITS) | mant;
r = uint32_as_float (res);
}
return r;
}
#define NBR_RND_MODES (4)
#define NBR_DAZ_MODES (2)
#define NBR_FTZ_MODES (2)
int main (void)
{
const uint32_t rnd_mode [NBR_RND_MODES] =
{
_MM_ROUND_NEAREST,
_MM_ROUND_TOWARD_ZERO,
_MM_ROUND_DOWN,
_MM_ROUND_UP
};
const uint32_t ftz_mode [NBR_FTZ_MODES] =
{
_MM_FLUSH_ZERO_OFF,
_MM_FLUSH_ZERO_ON
};
const uint32_t daz_mode [NBR_DAZ_MODES] =
{
_MM_DENORMALS_ZERO_OFF,
_MM_DENORMALS_ZERO_ON
};
uint32_t iarg, ires, iref;
float arg, res, ref;
double relerr, maxrelerr;
for (int rnd = 0; rnd < NBR_RND_MODES; rnd++) {
printf ("rnd=%d\n", rnd);
for (int ftz = 0; ftz < NBR_FTZ_MODES; ftz++) {
printf (" ftz=%d\n", ftz);
for (int daz = 0; daz < NBR_DAZ_MODES; daz++) {
printf (" daz=%d\n", daz); fflush(stdout);
maxrelerr = 0;
iarg = 0;
do {
arg = uint32_as_float (iarg);
ref = sse_rsqrtf (arg, daz_mode[daz], ftz_mode[ftz], rnd_mode[rnd]);
res = my_rsqrtf (arg);
if ((arg >= 1.17549435e-38f) && (arg < 3.40282347e+38f)) { /* normals only */
relerr = fabs ((ref - sqrt(1.0/(double)arg)) / sqrt(1.0/(double)arg));
if (relerr > maxrelerr) maxrelerr = relerr;
}
iref = float_as_uint32 (ref);
ires = float_as_uint32 (res);
if (ires != iref) {
printf ("!!!! rnd=%d ftz=%d daz=%d arg=%08x res=%08x ref=%08x\n",
rnd, ftz, daz, iarg, ires, iref);
return EXIT_FAILURE;
}
iarg++;
} while (iarg);
printf (" maxrelerr = %15.8e\n", maxrelerr);
}
}
}
printf ("RSQRTSS emulation test passed\n");
return EXIT_SUCCESS;
}
我正在寻找一个 C++ 函数，它 returns 一个浮点数的逆平方根：rsqrt(x) = 1/sqrt(x) 通过使用像内置 XMM 操作 RSQRTSS 这样的精确方法（cf. https://www.felixcloutier.com/x86/rsqrtss)。 （即，我想要内置的近似值而不是更精确的 1/sqrtf，而且我不关心速度（很多）。）
根据这个问题：
Is there a fast C or C++ standard library function for double precision inverse square root?
...至少没有 "fast way with double precision" 可以使用标准 C++ 库来完成此操作。但是如何做到慢，非标准和 floats？
标准库中没有执行此操作的函数，但您的编译器可能会优化表达式 1 / sqrt(value)，使其发出 RSQRTSS 指令。例如，使用编译器标志 -ffast-math -march=native，GCC 将发出该指令，请参阅：https://godbolt.org/z/cL6seG
就其价值而言，正如@François Andrieux 所建议的那样，我最终在 C++ 中以普通汇编实现了它（更准确地说，我使用了 ASMJIT）。
这很好用，尽管它有失去可移植性的缺点（虽然不如普通 asm）。但这在某种程度上是我的问题所固有的，因为我想使用一个非常具体的 x86 函数。
这是我的代码：
typedef float(*JITFunc)();
JITFunc func;
asmjit::JitRuntime jit_runtime;
asmjit::CodeHolder code;
code.init(jit_runtime.getCodeInfo());
asmjit::X86Compiler cc(&code);
cc.addFunc(asmjit::FuncSignature0<float>());
float value = 2.71; // Some example value.
asmjit::X86Xmm x = cc.newXmm();
setXmmVar(cc, x, value);
cc.rsqrtss(x, x); // THE asm function.
cc.ret(x);
cc.endFunc();
cc.finalize();
jit_runtime.add(&func, &code);
return func(); // Or something to that effect. func() is the result, anyway.
RSQRTSS 指令可通过 immintrin.h 中声明的 _mm_rsqrt_ss() 内在函数轻松访问。但是我们也可以在软件中模拟指令，就像下面的 my_rsqrtf() 函数所做的那样。通过简单地观察 RSQRTSS 的输出，很容易发现它的函数值基于 211 个条目的（虚拟）table，每个条目的大小为 12 位.
注意 "virtual" 属性，因为硬件不太可能使用直接 24 Kbit table。我对 table 条目中模式的分析不建议使用 bipartite tables。一种更简单的压缩方案——就像我在下面的代码中使用的那样——基于 table 基值和 table 偏移量 可能 是用过的。我的方案只需要一个窄加器，但将 ROM 存储减少到 13 Kbit，即几乎一半。
下面的实现是在使用 Ivy Bridge 架构的英特尔至强处理器 E3-1270 V2 上开发和测试的。 RSQRTSS 在各种 Intel 架构之间的实现可能存在一些功能差异，并且这种差异可能存在于来自不同 x86-84 供应商的架构之间。
下面的框架检查 my_rsqrtf() 的仿真为所有四种舍入模式、两种 DAZ（非正规数为零）模式和两种 FTZ（齐平）提供与 RSQRTSS 相同的结果到零）模式。我们发现函数结果不受任何模式的影响，这与英特尔在 RSQRTSS 中指定的方式相符
英特尔® 64 位和 IA-32 架构软件开发人员手册：
The RSQRTSS instruction is not affected by the rounding control bits in the MXCSR register. When a source value is a 0.0, an ∞ of the sign of the source value is returned. A denormal source value is treated as a 0.0 (of the same sign). When a source value is a negative value (other than −0.0), a floating-point indefinite is returned. When a source value is an SNaN or QNaN, the SNaN is converted to a QNaN or the source QNaN is returned.
#include <stdio.h>
#include <stdlib.h>
#include <stdint.h>
#include <string.h>
#include <math.h>
#include "immintrin.h"
/* SSE reference for RSQRTSS instruction */
float sse_rsqrtf (float a, uint32_t daz, uint32_t ftz, uint32_t rnd)
{
__m128 b, t;
float res;
uint32_t old_mxcsr;
old_mxcsr = _mm_getcsr();
_MM_SET_DENORMALS_ZERO_MODE (daz);
_MM_SET_FLUSH_ZERO_MODE (ftz);
_MM_SET_ROUNDING_MODE (rnd);
b = _mm_set_ss (a);
t = _mm_rsqrt_ss (b);
_mm_store_ss (&res, t);
_mm_setcsr (old_mxcsr);
return res;
}
inline uint32_t float_as_uint32 (float a)
{
uint32_t r;
memcpy (&r, &a, sizeof r);
return r;
}
inline float uint32_as_float (uint32_t a)
{
float r;
memcpy (&r, &a, sizeof r);
return r;
}
#define LOG2_NBR_TAB_ENTRIES (11)
#define NBR_TAB_ENTRIES (1 << LOG2_NBR_TAB_ENTRIES)
#define TAB_ENTRY_BITS (12)
#define BASE_TAB_ENTRY_BITS (9)
/* 128 9-bit entries = 1152 bits */
const uint16_t base_tab[128] = {
0x0ce, 0x0c9, 0x0c3, 0x0be, 0x0b9, 0x0b4, 0x0af, 0x0aa,
0x0a6, 0x0a1, 0x09d, 0x098, 0x094, 0x090, 0x08b, 0x087,
0x083, 0x07f, 0x07c, 0x078, 0x074, 0x070, 0x06d, 0x069,
0x066, 0x062, 0x05f, 0x05c, 0x058, 0x055, 0x052, 0x04f,
0x04c, 0x049, 0x046, 0x043, 0x040, 0x03d, 0x03a, 0x038,
0x035, 0x032, 0x030, 0x02d, 0x02a, 0x028, 0x025, 0x023,
0x021, 0x01e, 0x01c, 0x019, 0x017, 0x015, 0x013, 0x010,
0x00e, 0x00c, 0x00a, 0x008, 0x006, 0x004, 0x002, 0x000,
0x1f8, 0x1f0, 0x1e9, 0x1e1, 0x1da, 0x1d3, 0x1cc, 0x1c5,
0x1bf, 0x1b8, 0x1b2, 0x1ab, 0x1a5, 0x19f, 0x199, 0x194,
0x18e, 0x188, 0x183, 0x17e, 0x178, 0x173, 0x16e, 0x169,
0x164, 0x15f, 0x15a, 0x156, 0x151, 0x14d, 0x148, 0x144,
0x13f, 0x13b, 0x137, 0x133, 0x12f, 0x12b, 0x127, 0x123,
0x11f, 0x11b, 0x118, 0x114, 0x110, 0x10d, 0x109, 0x106,
0x102, 0x0ff, 0x0fc, 0x0f8, 0x0f5, 0x0f2, 0x0ef, 0x0eb,
0x0e8, 0x0e5, 0x0e2, 0x0df, 0x0dc, 0x0d9, 0x0d7, 0x0d4,
};
/* 2048 6-bit entries = 12288 bits */
const uint8_t ofs_tab[2048] = {
0x2f, 0x2c, 0x2a, 0x27, 0x24, 0x21, 0x1e, 0x1c, 0x19, 0x16, 0x13, 0x10, 0x0e, 0x0b, 0x08, 0x05,
0x2b, 0x28, 0x25, 0x22, 0x1f, 0x1d, 0x1a, 0x17, 0x15, 0x12, 0x0f, 0x0c, 0x0a, 0x07, 0x04, 0x02,
0x2f, 0x2c, 0x29, 0x27, 0x24, 0x21, 0x1f, 0x1c, 0x19, 0x17, 0x14, 0x11, 0x0f, 0x0c, 0x09, 0x07,
0x2c, 0x29, 0x27, 0x24, 0x22, 0x1f, 0x1c, 0x1a, 0x17, 0x15, 0x12, 0x0f, 0x0d, 0x0a, 0x08, 0x05,
0x2a, 0x28, 0x25, 0x23, 0x20, 0x1e, 0x1b, 0x18, 0x16, 0x13, 0x11, 0x0e, 0x0c, 0x09, 0x07, 0x04,
0x2a, 0x27, 0x24, 0x22, 0x1f, 0x1d, 0x1a, 0x18, 0x15, 0x13, 0x10, 0x0e, 0x0b, 0x09, 0x07, 0x04,
0x2a, 0x27, 0x25, 0x22, 0x20, 0x1d, 0x1b, 0x18, 0x16, 0x13, 0x11, 0x0f, 0x0c, 0x0a, 0x07, 0x05,
0x2a, 0x28, 0x26, 0x23, 0x21, 0x1e, 0x1c, 0x1a, 0x17, 0x15, 0x12, 0x10, 0x0e, 0x0b, 0x09, 0x07,
0x24, 0x22, 0x1f, 0x1d, 0x1b, 0x18, 0x16, 0x14, 0x11, 0x0f, 0x0d, 0x0a, 0x08, 0x06, 0x03, 0x01,
0x27, 0x24, 0x22, 0x20, 0x1d, 0x1b, 0x19, 0x16, 0x14, 0x12, 0x10, 0x0d, 0x0b, 0x09, 0x06, 0x04,
0x22, 0x20, 0x1d, 0x1b, 0x19, 0x17, 0x14, 0x12, 0x10, 0x0e, 0x0b, 0x09, 0x07, 0x05, 0x02, 0x00,
0x26, 0x24, 0x21, 0x1f, 0x1d, 0x1b, 0x19, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0b, 0x09, 0x07, 0x05,
0x23, 0x20, 0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x13, 0x11, 0x0f, 0x0d, 0x0b, 0x09, 0x06, 0x04, 0x02,
0x20, 0x1e, 0x1c, 0x1a, 0x17, 0x15, 0x13, 0x11, 0x0f, 0x0d, 0x0b, 0x09, 0x06, 0x04, 0x02, 0x00,
0x26, 0x24, 0x22, 0x20, 0x1e, 0x1c, 0x19, 0x17, 0x15, 0x13, 0x11, 0x0f, 0x0d, 0x0b, 0x09, 0x07,
0x25, 0x23, 0x21, 0x1f, 0x1d, 0x1a, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x08, 0x06,
0x24, 0x22, 0x20, 0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x08, 0x06,
0x24, 0x22, 0x20, 0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x08, 0x06,
0x1d, 0x1b, 0x19, 0x17, 0x15, 0x13, 0x11, 0x0f, 0x0d, 0x0b, 0x09, 0x07, 0x05, 0x03, 0x01, 0x00,
0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0b, 0x09, 0x07, 0x05, 0x03, 0x01,
0x1f, 0x1d, 0x1c, 0x1a, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0d, 0x0b, 0x09, 0x07, 0x05, 0x03,
0x21, 0x20, 0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x15, 0x13, 0x11, 0x0f, 0x0d, 0x0b, 0x0a, 0x08, 0x06,
0x1c, 0x1a, 0x19, 0x17, 0x15, 0x13, 0x11, 0x10, 0x0e, 0x0c, 0x0a, 0x08, 0x07, 0x05, 0x03, 0x01,
0x1f, 0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x15, 0x13, 0x11, 0x0f, 0x0e, 0x0c, 0x0a, 0x08, 0x07, 0x05,
0x1b, 0x19, 0x18, 0x16, 0x14, 0x12, 0x11, 0x0f, 0x0d, 0x0b, 0x0a, 0x08, 0x06, 0x04, 0x03, 0x01,
0x1f, 0x1e, 0x1c, 0x1a, 0x18, 0x17, 0x15, 0x13, 0x12, 0x10, 0x0e, 0x0c, 0x0b, 0x09, 0x07, 0x06,
0x1c, 0x1a, 0x19, 0x17, 0x15, 0x13, 0x12, 0x10, 0x0e, 0x0d, 0x0b, 0x09, 0x08, 0x06, 0x04, 0x03,
0x19, 0x17, 0x16, 0x14, 0x12, 0x11, 0x0f, 0x0d, 0x0c, 0x0a, 0x08, 0x07, 0x05, 0x03, 0x02, 0x00,
0x1f, 0x1d, 0x1b, 0x1a, 0x18, 0x16, 0x15, 0x13, 0x11, 0x10, 0x0e, 0x0d, 0x0b, 0x09, 0x08, 0x06,
0x1d, 0x1b, 0x19, 0x18, 0x16, 0x14, 0x13, 0x11, 0x10, 0x0e, 0x0c, 0x0b, 0x09, 0x08, 0x06, 0x04,
0x1b, 0x19, 0x18, 0x16, 0x15, 0x13, 0x11, 0x10, 0x0e, 0x0d, 0x0b, 0x0a, 0x08, 0x06, 0x05, 0x03,
0x1a, 0x18, 0x17, 0x15, 0x13, 0x12, 0x10, 0x0f, 0x0d, 0x0c, 0x0a, 0x09, 0x07, 0x06, 0x04, 0x02,
0x19, 0x17, 0x16, 0x14, 0x13, 0x11, 0x10, 0x0e, 0x0d, 0x0b, 0x0a, 0x08, 0x07, 0x05, 0x03, 0x02,
0x18, 0x17, 0x15, 0x14, 0x12, 0x11, 0x0f, 0x0e, 0x0c, 0x0b, 0x09, 0x08, 0x06, 0x05, 0x03, 0x02,
0x18, 0x17, 0x15, 0x14, 0x12, 0x11, 0x0f, 0x0e, 0x0d, 0x0b, 0x0a, 0x08, 0x07, 0x05, 0x04, 0x02,
0x19, 0x17, 0x16, 0x14, 0x13, 0x11, 0x10, 0x0e, 0x0d, 0x0c, 0x0a, 0x09, 0x07, 0x06, 0x04, 0x03,
0x19, 0x18, 0x16, 0x15, 0x14, 0x12, 0x11, 0x0f, 0x0e, 0x0c, 0x0b, 0x0a, 0x08, 0x07, 0x05, 0x04,
0x1a, 0x19, 0x18, 0x16, 0x15, 0x13, 0x12, 0x10, 0x0f, 0x0e, 0x0c, 0x0b, 0x09, 0x08, 0x07, 0x05,
0x1c, 0x1a, 0x19, 0x18, 0x16, 0x15, 0x13, 0x12, 0x11, 0x0f, 0x0e, 0x0c, 0x0b, 0x0a, 0x08, 0x07,
0x15, 0x14, 0x13, 0x11, 0x10, 0x0f, 0x0d, 0x0c, 0x0a, 0x09, 0x08, 0x06, 0x05, 0x04, 0x02, 0x01,
0x17, 0x16, 0x15, 0x13, 0x12, 0x11, 0x0f, 0x0e, 0x0d, 0x0b, 0x0a, 0x08, 0x07, 0x06, 0x04, 0x03,
0x1a, 0x18, 0x17, 0x16, 0x14, 0x13, 0x12, 0x10, 0x0f, 0x0e, 0x0c, 0x0b, 0x0a, 0x08, 0x07, 0x06,
0x14, 0x13, 0x12, 0x10, 0x0f, 0x0e, 0x0c, 0x0b, 0x0a, 0x08, 0x07, 0x06, 0x05, 0x03, 0x02, 0x01,
0x17, 0x16, 0x15, 0x13, 0x12, 0x11, 0x0f, 0x0e, 0x0d, 0x0c, 0x0a, 0x09, 0x08, 0x06, 0x05, 0x04,
0x1b, 0x19, 0x18, 0x17, 0x15, 0x14, 0x13, 0x12, 0x10, 0x0f, 0x0e, 0x0c, 0x0b, 0x0a, 0x09, 0x07,
0x16, 0x15, 0x13, 0x12, 0x11, 0x10, 0x0e, 0x0d, 0x0c, 0x0b, 0x09, 0x08, 0x07, 0x06, 0x04, 0x03,
0x1a, 0x19, 0x17, 0x16, 0x15, 0x14, 0x12, 0x11, 0x10, 0x0f, 0x0d, 0x0c, 0x0b, 0x0a, 0x08, 0x07,
0x16, 0x15, 0x13, 0x12, 0x11, 0x10, 0x0e, 0x0d, 0x0c, 0x0b, 0x0a, 0x08, 0x07, 0x06, 0x05, 0x03,
0x12, 0x11, 0x10, 0x0f, 0x0d, 0x0c, 0x0b, 0x0a, 0x08, 0x07, 0x06, 0x05, 0x04, 0x02, 0x01, 0x00,
0x17, 0x16, 0x14, 0x13, 0x12, 0x11, 0x10, 0x0e, 0x0d, 0x0c, 0x0b, 0x0a, 0x08, 0x07, 0x06, 0x05,
0x14, 0x12, 0x11, 0x10, 0x0f, 0x0e, 0x0d, 0x0b, 0x0a, 0x09, 0x08, 0x07, 0x05, 0x04, 0x03, 0x02,
0x19, 0x18, 0x16, 0x15, 0x14, 0x13, 0x12, 0x11, 0x0f, 0x0e, 0x0d, 0x0c, 0x0b, 0x0a, 0x08, 0x07,
0x16, 0x15, 0x14, 0x13, 0x11, 0x10, 0x0f, 0x0e, 0x0d, 0x0c, 0x0b, 0x09, 0x08, 0x07, 0x06, 0x05,
0x14, 0x13, 0x11, 0x10, 0x0f, 0x0e, 0x0d, 0x0c, 0x0b, 0x09, 0x08, 0x07, 0x06, 0x05, 0x04, 0x03,
0x11, 0x10, 0x0f, 0x0e, 0x0d, 0x0c, 0x0b, 0x0a, 0x08, 0x07, 0x06, 0x05, 0x04, 0x03, 0x02, 0x01,
0x18, 0x16, 0x15, 0x14, 0x13, 0x12, 0x11, 0x10, 0x0f, 0x0e, 0x0c, 0x0b, 0x0a, 0x09, 0x08, 0x07,
0x16, 0x15, 0x14, 0x12, 0x11, 0x10, 0x0f, 0x0e, 0x0d, 0x0c, 0x0b, 0x0a, 0x09, 0x08, 0x06, 0x05,
0x14, 0x13, 0x12, 0x11, 0x10, 0x0f, 0x0e, 0x0d, 0x0c, 0x0b, 0x09, 0x08, 0x07, 0x06, 0x05, 0x04,
0x13, 0x12, 0x11, 0x10, 0x0f, 0x0e, 0x0d, 0x0b, 0x0a, 0x09, 0x08, 0x07, 0x06, 0x05, 0x04, 0x03,
0x12, 0x11, 0x10, 0x0f, 0x0e, 0x0d, 0x0c, 0x0a, 0x09, 0x08, 0x07, 0x06, 0x05, 0x04, 0x03, 0x02,
0x11, 0x10, 0x0f, 0x0e, 0x0d, 0x0c, 0x0b, 0x0a, 0x09, 0x08, 0x07, 0x06, 0x04, 0x03, 0x02, 0x01,
0x10, 0x0f, 0x0e, 0x0d, 0x0c, 0x0b, 0x0a, 0x09, 0x08, 0x07, 0x06, 0x05, 0x04, 0x03, 0x02, 0x01,
0x10, 0x0f, 0x0e, 0x0d, 0x0c, 0x0b, 0x0a, 0x09, 0x08, 0x07, 0x06, 0x05, 0x04, 0x03, 0x02, 0x01,
0x10, 0x0f, 0x0e, 0x0d, 0x0c, 0x0b, 0x0a, 0x09, 0x08, 0x07, 0x06, 0x05, 0x04, 0x03, 0x02, 0x01,
0x3e, 0x3a, 0x36, 0x32, 0x2e, 0x2a, 0x26, 0x22, 0x1e, 0x1a, 0x16, 0x12, 0x0e, 0x0b, 0x07, 0x03,
0x3f, 0x3b, 0x37, 0x33, 0x2f, 0x2b, 0x27, 0x24, 0x20, 0x1c, 0x18, 0x14, 0x10, 0x0c, 0x09, 0x05,
0x39, 0x35, 0x31, 0x2e, 0x2a, 0x26, 0x22, 0x1e, 0x1b, 0x17, 0x13, 0x0f, 0x0c, 0x08, 0x04, 0x00,
0x3d, 0x39, 0x35, 0x31, 0x2e, 0x2a, 0x26, 0x23, 0x1f, 0x1b, 0x18, 0x14, 0x10, 0x0d, 0x09, 0x05,
0x3a, 0x36, 0x32, 0x2f, 0x2b, 0x27, 0x24, 0x20, 0x1d, 0x19, 0x15, 0x12, 0x0e, 0x0b, 0x07, 0x03,
0x38, 0x34, 0x31, 0x2d, 0x2a, 0x26, 0x22, 0x1f, 0x1b, 0x18, 0x14, 0x11, 0x0d, 0x0a, 0x06, 0x03,
0x37, 0x34, 0x30, 0x2d, 0x29, 0x26, 0x22, 0x1f, 0x1b, 0x18, 0x15, 0x11, 0x0e, 0x0a, 0x07, 0x03,
0x38, 0x35, 0x31, 0x2e, 0x2a, 0x27, 0x24, 0x20, 0x1d, 0x19, 0x16, 0x13, 0x0f, 0x0c, 0x09, 0x05,
0x32, 0x2e, 0x2b, 0x28, 0x24, 0x21, 0x1e, 0x1a, 0x17, 0x14, 0x11, 0x0d, 0x0a, 0x07, 0x03, 0x00,
0x35, 0x31, 0x2e, 0x2b, 0x28, 0x24, 0x21, 0x1e, 0x1b, 0x17, 0x14, 0x11, 0x0e, 0x0a, 0x07, 0x04,
0x31, 0x2e, 0x2a, 0x27, 0x24, 0x21, 0x1e, 0x1a, 0x17, 0x14, 0x11, 0x0e, 0x0b, 0x07, 0x04, 0x01,
0x36, 0x33, 0x30, 0x2c, 0x29, 0x26, 0x23, 0x20, 0x1d, 0x1a, 0x17, 0x13, 0x10, 0x0d, 0x0a, 0x07,
0x34, 0x31, 0x2e, 0x2b, 0x28, 0x25, 0x21, 0x1e, 0x1b, 0x18, 0x15, 0x12, 0x0f, 0x0c, 0x09, 0x06,
0x33, 0x30, 0x2d, 0x2a, 0x27, 0x24, 0x21, 0x1e, 0x1b, 0x18, 0x15, 0x12, 0x0f, 0x0c, 0x09, 0x06,
0x33, 0x30, 0x2d, 0x2a, 0x27, 0x24, 0x21, 0x1e, 0x1b, 0x18, 0x15, 0x13, 0x10, 0x0d, 0x0a, 0x07,
0x2c, 0x29, 0x26, 0x23, 0x20, 0x1d, 0x1a, 0x18, 0x15, 0x12, 0x0f, 0x0c, 0x09, 0x06, 0x03, 0x01,
0x2e, 0x2b, 0x28, 0x25, 0x22, 0x1f, 0x1d, 0x1a, 0x17, 0x14, 0x11, 0x0e, 0x0c, 0x09, 0x06, 0x03,
0x30, 0x2e, 0x2b, 0x28, 0x25, 0x22, 0x20, 0x1d, 0x1a, 0x17, 0x14, 0x12, 0x0f, 0x0c, 0x09, 0x07,
0x2c, 0x29, 0x26, 0x24, 0x21, 0x1e, 0x1b, 0x19, 0x16, 0x13, 0x10, 0x0e, 0x0b, 0x08, 0x06, 0x03,
0x28, 0x25, 0x23, 0x20, 0x1d, 0x1b, 0x18, 0x15, 0x13, 0x10, 0x0d, 0x0b, 0x08, 0x05, 0x03, 0x00,
0x2d, 0x2b, 0x28, 0x25, 0x23, 0x20, 0x1d, 0x1b, 0x18, 0x15, 0x13, 0x10, 0x0e, 0x0b, 0x08, 0x06,
0x2b, 0x28, 0x26, 0x23, 0x21, 0x1e, 0x1b, 0x19, 0x16, 0x14, 0x11, 0x0f, 0x0c, 0x09, 0x07, 0x04,
0x2a, 0x27, 0x25, 0x22, 0x1f, 0x1d, 0x1a, 0x18, 0x15, 0x13, 0x10, 0x0e, 0x0b, 0x09, 0x06, 0x03,
0x29, 0x26, 0x24, 0x21, 0x1f, 0x1c, 0x1a, 0x17, 0x15, 0x12, 0x10, 0x0d, 0x0b, 0x08, 0x06, 0x03,
0x29, 0x26, 0x24, 0x21, 0x1f, 0x1d, 0x1a, 0x18, 0x15, 0x13, 0x10, 0x0e, 0x0b, 0x09, 0x06, 0x04,
0x2a, 0x27, 0x25, 0x22, 0x20, 0x1d, 0x1b, 0x19, 0x16, 0x14, 0x11, 0x0f, 0x0d, 0x0a, 0x08, 0x05,
0x2b, 0x29, 0x26, 0x24, 0x21, 0x1f, 0x1d, 0x1a, 0x18, 0x15, 0x13, 0x11, 0x0e, 0x0c, 0x0a, 0x07,
0x25, 0x23, 0x20, 0x1e, 0x1b, 0x19, 0x17, 0x14, 0x12, 0x10, 0x0d, 0x0b, 0x09, 0x06, 0x04, 0x02,
0x27, 0x25, 0x23, 0x20, 0x1e, 0x1c, 0x1a, 0x17, 0x15, 0x13, 0x10, 0x0e, 0x0c, 0x09, 0x07, 0x05,
0x23, 0x20, 0x1e, 0x1c, 0x1a, 0x17, 0x15, 0x13, 0x10, 0x0e, 0x0c, 0x0a, 0x07, 0x05, 0x03, 0x01,
0x26, 0x24, 0x22, 0x20, 0x1d, 0x1b, 0x19, 0x17, 0x15, 0x12, 0x10, 0x0e, 0x0c, 0x09, 0x07, 0x05,
0x23, 0x21, 0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x13, 0x11, 0x0f, 0x0d, 0x0b, 0x08, 0x06, 0x04, 0x02,
0x28, 0x25, 0x23, 0x21, 0x1f, 0x1d, 0x1b, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x09, 0x07,
0x25, 0x23, 0x21, 0x1f, 0x1d, 0x1a, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x09, 0x07, 0x05,
0x23, 0x21, 0x1f, 0x1d, 0x1b, 0x19, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x08, 0x06, 0x04,
0x22, 0x1f, 0x1d, 0x1b, 0x19, 0x17, 0x15, 0x13, 0x11, 0x0f, 0x0d, 0x0b, 0x09, 0x07, 0x05, 0x03,
0x21, 0x1f, 0x1c, 0x1a, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x08, 0x06, 0x04, 0x02,
0x20, 0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x08, 0x06, 0x04, 0x02,
0x20, 0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x08, 0x06, 0x04, 0x02,
0x20, 0x1f, 0x1d, 0x1b, 0x19, 0x17, 0x15, 0x13, 0x11, 0x0f, 0x0d, 0x0b, 0x09, 0x07, 0x05, 0x03,
0x21, 0x1f, 0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x08, 0x07, 0x05,
0x23, 0x21, 0x1f, 0x1d, 0x1b, 0x19, 0x17, 0x15, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x08, 0x06,
0x1c, 0x1b, 0x19, 0x17, 0x15, 0x13, 0x11, 0x0f, 0x0e, 0x0c, 0x0a, 0x08, 0x06, 0x04, 0x02, 0x01,
0x1f, 0x1d, 0x1b, 0x19, 0x17, 0x15, 0x14, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x09, 0x07, 0x05, 0x03,
0x21, 0x1f, 0x1e, 0x1c, 0x1a, 0x18, 0x16, 0x15, 0x13, 0x11, 0x0f, 0x0d, 0x0c, 0x0a, 0x08, 0x06,
0x1c, 0x1b, 0x19, 0x17, 0x15, 0x13, 0x12, 0x10, 0x0e, 0x0c, 0x0a, 0x09, 0x07, 0x05, 0x03, 0x02,
0x20, 0x1e, 0x1c, 0x1a, 0x19, 0x17, 0x15, 0x13, 0x12, 0x10, 0x0e, 0x0c, 0x0b, 0x09, 0x07, 0x05,
0x1c, 0x1a, 0x18, 0x16, 0x15, 0x13, 0x11, 0x0f, 0x0e, 0x0c, 0x0a, 0x08, 0x07, 0x05, 0x03, 0x01,
0x20, 0x1e, 0x1c, 0x1b, 0x19, 0x17, 0x15, 0x14, 0x12, 0x10, 0x0f, 0x0d, 0x0b, 0x09, 0x08, 0x06,
0x1c, 0x1b, 0x19, 0x17, 0x15, 0x14, 0x12, 0x10, 0x0f, 0x0d, 0x0b, 0x0a, 0x08, 0x06, 0x05, 0x03,
0x19, 0x17, 0x16, 0x14, 0x12, 0x11, 0x0f, 0x0d, 0x0c, 0x0a, 0x08, 0x07, 0x05, 0x03, 0x02, 0x00,
0x1e, 0x1d, 0x1b, 0x19, 0x18, 0x16, 0x14, 0x13, 0x11, 0x10, 0x0e, 0x0c, 0x0b, 0x09, 0x07, 0x06,
0x1c, 0x1a, 0x19, 0x17, 0x16, 0x14, 0x12, 0x11, 0x0f, 0x0d, 0x0c, 0x0a, 0x08, 0x07, 0x05, 0x04,
0x1a, 0x18, 0x17, 0x15, 0x14, 0x12, 0x10, 0x0f, 0x0d, 0x0b, 0x0a, 0x08, 0x07, 0x05, 0x03, 0x02,
0x18, 0x17, 0x15, 0x13, 0x12, 0x10, 0x0f, 0x0d, 0x0c, 0x0a, 0x08, 0x07, 0x05, 0x04, 0x02, 0x00,
0x1f, 0x1d, 0x1c, 0x1a, 0x19, 0x17, 0x15, 0x14, 0x12, 0x11, 0x0f, 0x0e, 0x0c, 0x0a, 0x09, 0x07,
0x1e, 0x1c, 0x1b, 0x19, 0x18, 0x16, 0x14, 0x13, 0x11, 0x10, 0x0e, 0x0d, 0x0b, 0x0a, 0x08, 0x07,
0x1d, 0x1c, 0x1a, 0x18, 0x17, 0x15, 0x14, 0x12, 0x11, 0x0f, 0x0e, 0x0c, 0x0b, 0x09, 0x08, 0x06,
0x1d, 0x1b, 0x1a, 0x18, 0x17, 0x15, 0x13, 0x12, 0x10, 0x0f, 0x0d, 0x0c, 0x0a, 0x09, 0x07, 0x06,
0x1c, 0x1b, 0x19, 0x18, 0x16, 0x15, 0x13, 0x12, 0x10, 0x0f, 0x0d, 0x0c, 0x0b, 0x09, 0x08, 0x06,
0x1d, 0x1b, 0x1a, 0x18, 0x17, 0x15, 0x14, 0x12, 0x11, 0x0f, 0x0e, 0x0c, 0x0b, 0x09, 0x08, 0x06,
0x1d, 0x1c, 0x1a, 0x19, 0x17, 0x16, 0x14, 0x13, 0x11, 0x10, 0x0e, 0x0d, 0x0b, 0x0a, 0x09, 0x07,
0x16, 0x14, 0x13, 0x11, 0x10, 0x0e, 0x0d, 0x0c, 0x0a, 0x09, 0x07, 0x06, 0x04, 0x03, 0x02, 0x00,
0x17, 0x15, 0x14, 0x12, 0x11, 0x10, 0x0e, 0x0d, 0x0b, 0x0a, 0x08, 0x07, 0x06, 0x04, 0x03, 0x01,
};
#define IEEE_BINARY32_EXPO_BIAS (127)
#define IEEE_BINARY32_MANT_BITS (23)
#define IEEE_BINARY32_EXPO_BITS (8)
#define IEEE_BINARY32_EXPO_MASK (0x7f800000)
#define IEEE_BINARY32_NAN_INDEF (0xffc00000)
#define IEEE_BINARY32_POS_INF (0x7f800000)
#define IEEE_BINARY32_POS_ZERO (0x00000000)
#define IEEE_BINARY32_MIN_NORM (0x00800000)
#define IEEE_BINARY32_SIGN_BIT (0x80000000)
/* Emulate the RSQRTSS instruction in software */
float my_rsqrtf (float x)
{
float r;
uint32_t arg, res, idx, expo, mant;
arg = float_as_uint32 (x);
/* zeros and subnormals */
if ((arg & ~IEEE_BINARY32_SIGN_BIT) < IEEE_BINARY32_MIN_NORM) {
res = IEEE_BINARY32_POS_INF | (arg & IEEE_BINARY32_SIGN_BIT);
r = uint32_as_float (res);
}
/* NaNs */
else if ((arg & ~IEEE_BINARY32_SIGN_BIT) > IEEE_BINARY32_POS_INF) {
r = x + x; // convert SNaN to QNaN
}
/* negative arguments */
else if (arg & IEEE_BINARY32_SIGN_BIT) {
res = IEEE_BINARY32_NAN_INDEF;
r = uint32_as_float (res);
}
/* positive infinity */
else if (arg == IEEE_BINARY32_POS_INF) {
res = IEEE_BINARY32_POS_ZERO;
r = uint32_as_float (res);
}
/* positive normals */
else {
/* extract exponent lsb and leading mantissa bits for table index */
expo = (arg & IEEE_BINARY32_EXPO_MASK) >> IEEE_BINARY32_MANT_BITS;
idx = (arg >> (IEEE_BINARY32_MANT_BITS - LOG2_NBR_TAB_ENTRIES + 1))
& (NBR_TAB_ENTRIES - 1);
/* compute exponent and mantissa of reciprocal square root */
expo = (3 * IEEE_BINARY32_EXPO_BIAS + ~expo) >> 1;
mant = (((base_tab [idx >> 4] << (TAB_ENTRY_BITS - BASE_TAB_ENTRY_BITS)) + ofs_tab [idx])
<< (IEEE_BINARY32_MANT_BITS - TAB_ENTRY_BITS));
/* combine exponent and mantissa bits to compute final result */
res = (expo << IEEE_BINARY32_MANT_BITS) | mant;
r = uint32_as_float (res);
}
return r;
}
#define NBR_RND_MODES (4)
#define NBR_DAZ_MODES (2)
#define NBR_FTZ_MODES (2)
int main (void)
{
const uint32_t rnd_mode [NBR_RND_MODES] =
{
_MM_ROUND_NEAREST,
_MM_ROUND_TOWARD_ZERO,
_MM_ROUND_DOWN,
_MM_ROUND_UP
};
const uint32_t ftz_mode [NBR_FTZ_MODES] =
{
_MM_FLUSH_ZERO_OFF,
_MM_FLUSH_ZERO_ON
};
const uint32_t daz_mode [NBR_DAZ_MODES] =
{
_MM_DENORMALS_ZERO_OFF,
_MM_DENORMALS_ZERO_ON
};
uint32_t iarg, ires, iref;
float arg, res, ref;
double relerr, maxrelerr;
for (int rnd = 0; rnd < NBR_RND_MODES; rnd++) {
printf ("rnd=%d\n", rnd);
for (int ftz = 0; ftz < NBR_FTZ_MODES; ftz++) {
printf (" ftz=%d\n", ftz);
for (int daz = 0; daz < NBR_DAZ_MODES; daz++) {
printf (" daz=%d\n", daz); fflush(stdout);
maxrelerr = 0;
iarg = 0;
do {
arg = uint32_as_float (iarg);
ref = sse_rsqrtf (arg, daz_mode[daz], ftz_mode[ftz], rnd_mode[rnd]);
res = my_rsqrtf (arg);
if ((arg >= 1.17549435e-38f) && (arg < 3.40282347e+38f)) { /* normals only */
relerr = fabs ((ref - sqrt(1.0/(double)arg)) / sqrt(1.0/(double)arg));
if (relerr > maxrelerr) maxrelerr = relerr;
}
iref = float_as_uint32 (ref);
ires = float_as_uint32 (res);
if (ires != iref) {
printf ("!!!! rnd=%d ftz=%d daz=%d arg=%08x res=%08x ref=%08x\n",
rnd, ftz, daz, iarg, ires, iref);
return EXIT_FAILURE;
}
iarg++;
} while (iarg);
printf (" maxrelerr = %15.8e\n", maxrelerr);
}
}
}
printf ("RSQRTSS emulation test passed\n");
return EXIT_SUCCESS;
}