如何在角度制下正确计算舍入的三角函数？

Question

如何定义以度为单位而非弧度为参数，并计算出这些参数的正确舍入结果的三角函数？

在将参数传递给弧度中相应的函数之前，将参数乘以M_PI/180.0并不起作用，因为M_PI/180.0不是π/180。《浮点运算手册》第5.5节提供了一种计算参数与π/180的正确舍入积的方法，但某些参数仍可能使得该积接近两个连续可表示浮点数的中间值，然后即使在弧度中应用正确舍入的函数也可能产生错误的最终结果。

两种可能有效的策略是使用更高精度和使用CRlibm中的sinpi、cospi、tanpi三角函数，它们分别计算sin(πx)、cos(πx)和tan(πx)。

对于后一种策略，仍然存在一个问题，即除以180不是对许多参数精确的。关于更高精度的策略（将参数乘以π/180的扩展精度表示，然后应用弧度的扩展精度函数），可能仍然存在“精确”的情况。该定理指出有理参数的sin、cos和tan的唯一有理结果仅在0中获得，这显然不适用于角度版本，如果对于某些浮点输入x，sindeg(x)恰好处于两个连续可表示浮点数的中点，则无论中间精度如何都不能保证最终结果正确舍入。

- Pascal Cuoq

Answer 1

唯一使得 cosdeg(360q) 有理数的有理数 q 的分母为1、2、3、4或6。Joerg Jahnel 的 这篇论文 在第六节中使用了域论的简短而美妙的证明。（事实上，作者使用欧拉函数对代数数 cosdeg(360q) 的次数进行了表征。）因此，不存在浮点数 q 使得 cosdeg(360q) 介于相邻两个浮点数之间。

所以我想答案是“大约和实现弧度的sin等函数的方式相同”，尽管 @gnasher729 提出了一个非常好的观点，即角度的参数减少要容易得多。

- tmyklebu

---

我会立即开始阅读论文，但是推断浮点数sindeg没有非平凡精确情况的结果不应该是“不存在浮点数q，使得sindeg（q）恰好处于两个相邻浮点数之间”吗？ - Pascal Cuoq

@PascalCuoq：是的。我（误）引用的结果表明，导致有理（不要紧邻两个浮点数之间）值的sindeg(q)或cosdeg(q)的唯一有理q是30度的倍数，答案总是无理数或1/2的倍数。这比您需要的结果更强，因为每个浮点数和每个中间浮点数都是有理数。 - tmyklebu

实际上，我更欣赏第三节中的初等证明（尽管它只涵盖了sindeg和cosdeg）。 - Pascal Cuoq

@PascalCuoq：如果你想得到类似的tan结果，http://oberlin.edu/faculty/jcalcut/tanpap.pdf看起来不错。 - Mark Dickinson

Answer 2

这很困难。但积极的一面是，你可以将参数减少到正负45度以内。因此，在正负45度之间需要正确舍入的结果。对于非常小的x，sin（x）约为x *（pi / 180），这已经足够难以精确舍入。

例如，要获得大部分正确舍入的正弦函数结果，请取-45 <= x <= 45。将x分成xhi = round（512 x）/ 512和xlo = x-xhi。让sin（x度）≈ax-bx^3。舍入a和b，使得计算s（x）a*xhi - b *（xhi^3）时精确计算。仔细计算余数sin（x度）- s（x）; 因为结果很小，所以舍入误差应该相当小。加上s（x），这将大多数时间给出正确舍入的结果。

- gnasher729

Answer 3

这是一个很难的问题。让我澄清一些要点：

• 需要什么精度的输出？是IEEE 754单精度、双精度还是非标准精度？此外，我假设输入（以度数表示的输入）应该和输出具有相同的精度，因为对于正常的弧度输入来说是这样的。
• 你的性能指标是什么？CRlibm经过优化，可以产生正确舍入的双精度结果。另一方面，MPFR用于任意精度，但当你只需要双精度输出时，它比CRlibm慢得多。
• 你的工作范围是什么？即[min argument，max argmunet]？对于CRlibm来说，这很重要，因为它适用于双精度范围。然而，对于MPFR来说，这并不真的重要。

如果你必须使用仅以度数为输入，请使用MPFR。让我提醒你，任何以度数表示的参数，当乘以（Pi/180）时，会产生一个超越数。然而，传递给三角函数的是浮点数表示，最好是四舍五入到最接近的整数，以工作精度。

我建议你按照以下步骤操作：

1. 尽可能使用C库，因为它比绑定库提供更好的性能。请使用MPFR。
2. 将MPFR精度设置得比目标精度高得多。例如（目标精度+300）。这样做可以避免任何精度损失，对于操作((Argument*Pi)/180)。在MPFR C库中可以轻松完成此操作，通过mpfr_set_default_prec()。
3. 进行操作：X_n=(Argument*Pi)/180，然后执行Sin(X_n)或者你想要的其他函数。在MPFR中有一个常量Pi，它在你的工作精度内表示。
4. 将结果舍入到目标精度。

由Muller编写的"Elementary functions"统计表明，如果工作精度略高于两倍目标精度，则大多数（而不是全部）难以处理的情况都能正确舍入。但是，在您的情况下，由于输入在理论上是超越的，为了安全起见，牺牲一些性能，使工作精度远高于目标精度。实际上，如果您需要双精度最终结果，10倍完全足够覆盖近乎100％的情况。

如果您需要低精度，即单精度或更低精度，则可以进行详尽的测试，以确定产生所有正确舍入的最低工作精度。

- Garp

---

如果存在可表示的浮点参数，其图像恰好处于两个可表示连续浮点值之间的中点位置，“将MPFR精度设置得比目标精度高得多”是行不通的。如果我在问题中没有表述清楚，对于“确切”情况可能仍然存在问题，我感到很抱歉。但是已经指出，对于sin、cos和tan，无需担心这种确切情况。 - Pascal Cuoq

2- “通过这样做，您可以避免操作（（参数* Pi）/ 180）的任何精度损失。” 由于MPFR是二进制的，会存在精度损失。例如，使用您的提案进行tandeg（90）计算将得到一个大但有限的值，并在目标精度中舍入为无穷大。正确的结果应该是NaN。 - Pascal Cuoq

是的，那是我的问题的一部分，tmyklebu已经回答了。然而，在他的回答之前，并不明显sindeg，cosdeg和tandeg没有任何确切的情况。请注意，我对sindeg，cosdeg和tandeg感兴趣。关于sin，cos，tan的定理并不直接适用。 - Pascal Cuoq

如果你想要特殊的输出正确无误，例如tan(90)，你应该手动处理它，因为tan(90)就像cos(90)或sin(0)一样是一个特殊情况。真正的90是一个超越数，当它传递给tan()或cos()时会产生特殊的输出。 - Garp

我不是在谈论sindeg，cosdeg等内容，我只是希望您能够以比目标更高的精度执行操作(deg*Pi)/180，以避免可能的精度损失。 - Garp

Answer 4

首先，您需要检测确切的情况，这已经有答案了。现在，对于其他情况，存在着著名的制表者困境。如果您的算术具有固定（且较小的）精度，并且希望得到中间精度所需的认证界限，那么已知有两种解决方案:

• 根据 Nesterenko 和 Waldschmidt 的定理获得一个基于界限，在我博士论文的第 4.3 节中有描述。（顺便说一句，我认为这也将为您提供确切情况的形式）。但是您会得到非常大的精度界限（至少几百万位？）。
• 找到最难处理的情况。仅需在[0,180]范围内搜索即可，因为任何更大的参数都将缩小到[0,180]范围内具有相同小数部分的值（因为周期是整数）。

- vinc17