fpu 中的小数字计算错误？
Fault in calculation with small numbers in fpu?
嗨，我正在尝试在 fpu 中执行此公式。
(y = v*t*sin(a) - 0.5*g*t^2)
我的 C++ 代码是：
typedef void(*Ipa_algorithm2)(double t, double alfa, double *return_value);
Ipa_algorithm2 count_y;
count_y = (Ipa_algorithm2)GetProcAddress(hInstLibrary, "ipa_algorithm");
t = t + 0.01; //going from cas = 0
(*count_y)(t,camera.angleY, &y); //t = cas;
我在 asm 中的代码是：
section .data
help_var dq 0
speed dq 40.0 ;v = rychlost
number dq 180.0
grav dq 4.906865 ;grav= 0,5*g
ipa_algorithm2:
push ebp
mov ebp, esp
finit
fld qword [speed]
fld qword [ebp+8]
fmul st1
fstp qword [help_var] ;v pomocny je v*t
fldpi
fld qword [ebp+16] ;na st0 je uhel a na st1 3,14
fmul st1 ;na st0 je uhel * 3,14
fld qword [number]
fxch st1 ;na st0 je uhel*3,14 na st1 je 180
fdiv st1 ;na st0 je uhel v radianech
fsin
fld qword [help_var]
fmul st1 ;na st0 je v*t*sin uhlu
fst qword [help_var]
finit
fld qword [ebp+8]
fld qword [ebp+8]
fmul st1
fld qword [grav]
fmul st1
fld qword [help_var]
fxch st1
fsub st1
mov eax,[ebp+24]
fstp qword [eax]
mov esp, ebp
pop ebp
ret 0
问题是，函数 ipa_algorithm2 从一开始就给我正确的数字（与在 C 中执行相同操作的程序的输出相比），但经过几个步骤后，结果开始变得越来越糟。我检查了 3 个小时的代码，没有发现任何错误。有没有可能是我数的数太小了，fpu 没法数了？
更新：，您在整个输入范围内得到的数字都是错误的，因此大概您在执行公式时遇到了一个常规错误，而不是特定于 FP 的舍入错误或数值 accuracy/stability 问题类型。对于给出错误答案的输入，在调试器中单步执行您的函数，然后查看寄存器值。
或者更好，用标量 AVX 指令重写它，因为标量 AVX 比 x87 和 更容易，所以工作标量实现是更好的起点。对于 sin()，调用矢量化 sin() 实现，或让 gcc 使用 -O3 -ffast-math 自动矢量化您的函数。
（参见 https://sourceware.org/glibc/wiki/libmvec：glibc 具有矢量化数学库函数。）
如果您最终想要快速运行的东西，那么从使用慢 fsin 指令的标量 x87 实现开始可能是最没有用的起点。对于您甚至不打算使用的指令集，良好干净的 C 语言比草率的 asm 实现要好。 （对于最终的优化版本，在大多数情况下，具有内在函数的 C 比手写的 asm 更有意义）。参见 http://agner.org/optimize/, and other links in the x86 tag wiki。
在游戏中存储角度
将方向存储为 [x,y] 矢量，而不是弧度角度。 （或学位）。使用归一化 xy 向量，将两个角度相加变成 2x2 矩阵乘法（通过旋转矩阵）。但是 sin 变得微不足道：如果你保持向量标准化 (x^2 + y^2 = 1.0) 那么 sin(angle) = angle.y.
尽可能避免使用实际角度，而是使用归一化向量。有时您需要 atan2，但通常很少，您可以只使用普通库版本。
如果您以数组结构格式存储您的 xy 对，它将对 SIMD 友好，并且您可以轻松地使用 8 个浮点 x 值和 8 个浮点匹配 y值。使用打包到单个 SIMD 向量中的方向向量来做事通常 NOT 最优；不要被 "vector".
这个词所迷惑
另请参阅 https://whosebug.com/tags/sse/info and especially SIMD at Insomniac Games (GDC 2015)。这将帮助您了解如何设计您的程序，以便您以后可以在值得的地方使用 SIMD 对其进行优化。 （您 没有 一开始就对所有内容进行矢量化，但更改数据布局通常需要大量工作，因此请首先考虑让您的数据 SIMD 友好。）
数值错误的可能来源（事实证明这不是真正的问题）
一个可能的原因：The worst-case error for the fsin instruction for small inputs is actually about 1.37 quintillion units in the last place, leaving fewer than four bits correct.。大多数现代数学库不使用 fsin 指令来计算 sin 函数，因为它速度不快并且某些输入的准确性很差。
此外，根据您构建代码的方式，有些东西（例如 MSVCRT 启动，如果您使用 Windows 并使用它的旧版本）may have set the x87 FPU to less than 80-bit precision（64 位尾数） .
你为什么用 asm 写这个？您需要有关如何提高效率的建议吗？您应该 return float in st0 作为 return 值，而不是通过指针 arg 存储。另外，不要使用 finit。我认为你这样做只是因为你在加载东西后没有平衡 x87 堆栈和 pops，所以在重复调用之后你会从 x87 堆栈溢出中获得 NaN。您仍然 return 在 returns void 的函数中使用 x87 堆栈非空，所以您仍然做错了并且可能会破坏调用者。
使用fstp 或fmulp 使堆栈保持平衡。使用 fld st0 而不是另一个负载。使用 fmul qword [grav_zrychleni] 而不是单独的 fld.
或者更好的是，使用 SSE2 或 AVX 进行标量双精度数学运算。除非你真的想要 80 位 long double.
嗨，我正在尝试在 fpu 中执行此公式。
(y = v*t*sin(a) - 0.5*g*t^2)
我的 C++ 代码是：
typedef void(*Ipa_algorithm2)(double t, double alfa, double *return_value);
Ipa_algorithm2 count_y;
count_y = (Ipa_algorithm2)GetProcAddress(hInstLibrary, "ipa_algorithm");
t = t + 0.01; //going from cas = 0
(*count_y)(t,camera.angleY, &y); //t = cas;
我在 asm 中的代码是：
section .data
help_var dq 0
speed dq 40.0 ;v = rychlost
number dq 180.0
grav dq 4.906865 ;grav= 0,5*g
ipa_algorithm2:
push ebp
mov ebp, esp
finit
fld qword [speed]
fld qword [ebp+8]
fmul st1
fstp qword [help_var] ;v pomocny je v*t
fldpi
fld qword [ebp+16] ;na st0 je uhel a na st1 3,14
fmul st1 ;na st0 je uhel * 3,14
fld qword [number]
fxch st1 ;na st0 je uhel*3,14 na st1 je 180
fdiv st1 ;na st0 je uhel v radianech
fsin
fld qword [help_var]
fmul st1 ;na st0 je v*t*sin uhlu
fst qword [help_var]
finit
fld qword [ebp+8]
fld qword [ebp+8]
fmul st1
fld qword [grav]
fmul st1
fld qword [help_var]
fxch st1
fsub st1
mov eax,[ebp+24]
fstp qword [eax]
mov esp, ebp
pop ebp
ret 0
问题是，函数 ipa_algorithm2 从一开始就给我正确的数字（与在 C 中执行相同操作的程序的输出相比），但经过几个步骤后，结果开始变得越来越糟。我检查了 3 个小时的代码，没有发现任何错误。有没有可能是我数的数太小了，fpu 没法数了？
更新：
或者更好，用标量 AVX 指令重写它，因为标量 AVX 比 x87 和 sin()，调用矢量化 sin() 实现，或让 gcc 使用 -O3 -ffast-math 自动矢量化您的函数。
（参见 https://sourceware.org/glibc/wiki/libmvec：glibc 具有矢量化数学库函数。）
如果您最终想要快速运行的东西，那么从使用慢 fsin 指令的标量 x87 实现开始可能是最没有用的起点。对于您甚至不打算使用的指令集，良好干净的 C 语言比草率的 asm 实现要好。 （对于最终的优化版本，在大多数情况下，具有内在函数的 C 比手写的 asm 更有意义）。参见 http://agner.org/optimize/, and other links in the x86 tag wiki。
在游戏中存储角度
将方向存储为 [x,y] 矢量，而不是弧度角度。 （或学位）。使用归一化 xy 向量，将两个角度相加变成 2x2 矩阵乘法（通过旋转矩阵）。但是 sin 变得微不足道：如果你保持向量标准化 (x^2 + y^2 = 1.0) 那么 sin(angle) = angle.y.
尽可能避免使用实际角度，而是使用归一化向量。有时您需要 atan2，但通常很少，您可以只使用普通库版本。
如果您以数组结构格式存储您的 xy 对，它将对 SIMD 友好，并且您可以轻松地使用 8 个浮点 x 值和 8 个浮点匹配 y值。使用打包到单个 SIMD 向量中的方向向量来做事通常 NOT 最优；不要被 "vector".
另请参阅 https://whosebug.com/tags/sse/info and especially SIMD at Insomniac Games (GDC 2015)。这将帮助您了解如何设计您的程序，以便您以后可以在值得的地方使用 SIMD 对其进行优化。 （您 没有 一开始就对所有内容进行矢量化，但更改数据布局通常需要大量工作，因此请首先考虑让您的数据 SIMD 友好。）
数值错误的可能来源（事实证明这不是真正的问题）
一个可能的原因：The worst-case error for the fsin instruction for small inputs is actually about 1.37 quintillion units in the last place, leaving fewer than four bits correct.。大多数现代数学库不使用 fsin 指令来计算 sin 函数，因为它速度不快并且某些输入的准确性很差。
此外，根据您构建代码的方式，有些东西（例如 MSVCRT 启动，如果您使用 Windows 并使用它的旧版本）may have set the x87 FPU to less than 80-bit precision（64 位尾数） .
你为什么用 asm 写这个？您需要有关如何提高效率的建议吗？您应该 return float in st0 作为 return 值，而不是通过指针 arg 存储。另外，不要使用 finit。我认为你这样做只是因为你在加载东西后没有平衡 x87 堆栈和 pops，所以在重复调用之后你会从 x87 堆栈溢出中获得 NaN。您仍然 return 在 returns void 的函数中使用 x87 堆栈非空，所以您仍然做错了并且可能会破坏调用者。
使用fstp 或fmulp 使堆栈保持平衡。使用 fld st0 而不是另一个负载。使用 fmul qword [grav_zrychleni] 而不是单独的 fld.
或者更好的是，使用 SSE2 或 AVX 进行标量双精度数学运算。除非你真的想要 80 位 long double.