问题标签 [numerical-stability]

我猜有一些我无法找到的标准技巧：无论如何，我想以数值稳定的方式计算一个非常接近 1 的数字的大幂（想想 1-p，其中 p<1e-17） . 1-p 在我的系统上被截断为 1。

使用对数的泰勒展开，我得到以下界限

我有什么更聪明的办法吗？

我有一个float x应该在 <0,1> 范围内的数字，但它经历了几次数值运算 - 结果可能稍微超出 <0,1>。

我需要将此结果转换为uint y使用整个范围的UInt32. 当然，我需要钳制x在 <0,1> 范围内并对其进行缩放。

但是哪种操作顺序更好？

或者

换句话说，最好先缩放，然后钳或钳，然后缩放？我对 IEEE 浮点表示不是很深，但我相信上述表达式的计算顺序是有区别的。

考虑以下代码：

实际上是这样的：

我也可以手动写成这样：

我的印象是乘法是一种更简单的运算，与除法运算相比，它不会受到舍入误差的影响。

另外，我有这种感觉，对于大多数日常人类数字，乘法运算会在除法之前产生一个“更大的数字”（假设使用的数字通常大于 1）。并且被划分后的较大数字在数值上更稳定。示例 .. 考虑5 * 7 / 2.3第一个操作 (mult) 的精确位置，因为这些数字以二进制精确表示。然后除法就完成了，它就像我们要得到的一样精确。但是考虑一下7 / 2.3 * 5，第一个操作是除法，我们已经产生了一个不能用二进制精确表示的数字，下一个操作（mult）通过乘法夸大了任何不精确性。

我的问题基本上是……这有关系吗？首先使用除法时我确实会失去精度，还是我完全安全地使用最适合我的任何操作顺序并且我会得到相同的结果？

我在一个大型 C 编程项目中遇到了一个数值问题。（这是统计研究，不是课堂作业）。第一步涉及计算 sqrt(x^2 + y) - x，我需要为正数，但有时即使 x > 0 且 y > 0，我也会得到 sqrt(x^2 + y) - x < 0。示例：

我的输出：

这个输出充满了奇怪的行为。强调：

1. 我分配了 y 0.0062529947，但它输出为 0.00625299476087093353271484375。
2. x*x + y 打印出与 x*x 相同的值。
3. sqrt(x*x + y) - x < 0。

为什么会发生1-3？

我应该提一下：我在具有 gcc 版本的 64 位 Mac OX 10.9.4 机器上运行了这个示例：

和具有 gcc 版本的 64 位 CentOS 服务器：

此外，编译在任何一台机器上都没有返回错误或警告：

在 GLSL（特别是我正在使用的 3.00）中，有两个版本 atan()：atan(y_over_x)只能返回 -PI/2、PI/2 之间的角度，而atan(y/x)可以考虑所有 4 个象限，因此角度范围涵盖了从 -PI 开始的所有内容， PI，很像atan2()C++。

我想使用第二个atan将 XY 坐标转换为角度。但是，atan()在 GLSL 中，除了不能处理 when 之外x = 0，还不是很稳定。尤其是在x接近零的情况下，除法可能会溢出，从而导致相反的结果角度（您会得到接近 -PI/2 的东西，而您假设会得到大约 PI/2）。

我们可以在 GLSL 之上构建一个好的、简单的实现atan(y,x)以使其更健壮吗？

我在 Matlab 中编程，在我的程序中我需要解决一个系统Ax=b，其中A是一个具有非常小的条目m的m方阵。如果m增加，则条目A变小。

A是一个稀疏矩阵，所以我用sparseMatlab的函数重命名这个矩阵：B=sparse(A).

我可以使用特殊的方法来解决Bx=b没有下溢或数字错误吗？

我正在用 C# 构建一个类似 CAD 的应用程序。我将SlimDX其用作图形引擎，对于数字运算部分，我构建了自定义库，这些库自然而然地最终依赖于 System.Math 类。

现在，问题是SlimDX库使用由浮点数据类型组成的结构，而 Math 类包含几个只接受和返回双精度对象的方法，例如：Math.Ceiling 和 Math.Sin。所以我发现自己不断地来回转换我的数据，从浮动到加倍。

这似乎不对。我并不关心演员阵容可能对性能造成的影响（也许我应该如此？），但担心他们可能导致的数值不稳定性，这更可怕。

所以我只是想知道您通常如何处理这些情况，因为我猜这一定不是罕见的情况。

对于熟悉该方法的人来说，众所周知，必须求解方程：

然而，F 通常不是线性的，因此得到的方程通常对 y(i+1) 有许多不同的解。我们正在寻找哪种解决方案，例如，对 Newton-Raphson 方法会做些什么，以便找到正确的零？任何帮助表示赞赏。

假设我们有一组浮点数，尾数为“m”，指数为“e”。假设我们想要逼近一个函数“f”。

从理论上我们知道，通常使用“范围缩减函数”，然后从这样的函数中我们得出全局函数值。

例如 let x = (sx,ex,mx) (sign exp and mantissa) then... log2(x) = ex + log2(1.mx) 所以基本上范围缩小函数是“log2(1.mx)” .

目前我已经实现了倒数、平方根、log2 和 exp2，最近我开始使用三角函数。但是，如果给定一个全局误差界限（尤其是 ulp 误差），我可能会为范围缩减函数推导出一个误差界限，是否有一些关于此类问题的研究？说到 log2(x)（例如），我想说...

“好的，我想要带有 k ulp 错误的 log2(x)，为了实现这一点，考虑到我们的浮点系统，我们需要用 p ulp 错误来近似 log2(1.mx)”

请记住，正如我所说，我们知道我们正在使用浮点数，但格式是通用的，所以它可能是经典的 F32，但即使例如 e=10，m = 8 end 等等。

我实际上找不到任何显示此类研究的参考资料。我拥有的参考资料（即穆勒书）没有以这种方式处理该主题，因此我正在寻找某种论文或类似的东西。你知道有什么参考吗？

我也试图自己推导出这样的约束，但这并不容易......

我了解数学上等效的算术运算如何由于数值错误（例如以不同顺序求和浮点数）而导致不同的结果。

然而，令我惊讶的是，添加零sum可以改变结果。我认为这始终适用于浮动，无论如何：x + 0. == x。

这是一个例子。我希望所有的行都完全为零。任何人都可以解释为什么会这样吗？

M和的较小值似乎不会发生Z。

我也确定了a.dtype==b.dtype。

这是另外一个示例，它还演示了 python 的内置sum函数按预期运行：

我正在使用 numpy V1.9.2。