【问题标题】:MatLab - variable precision arithmeticMatLab - 可变精度算术
【发布时间】:2012-06-03 18:47:18
【问题描述】:

我有一个关于 vpa 命令的简短问题,可用于在 MatLab 中评估符号表达式。

我的教科书是这样说的:

"在数字上使用诸如sqrt之类的函数时需要小心,默认情况下会生成双精度浮点数。您需要将此类输入作为符号字符串传递给vpa正确评价:vpa('sqrt(5)/pi')。”

我不太明白这里的行话。为什么对于大多数输入,无论我输入vpa(input) 还是vpa('input'),我都会得到完全相同的答案,但对于平方根却没有?例如,如果我输入vpa(sin(pi/4))vpa('sin(pi/4)'),我会得到完全相同的答案,但如果我输入上面给出的问题vpa(sqrt(5)/pi),我不会得到与输入vpa('sqrt(5)/pi') 相同的答案。

如果有人能比我上面的书更详细地解释这一点,我将不胜感激!

【问题讨论】:

标签: matlab numbers arbitrary-precision


【解决方案1】:

永远不要假设像 vpa(sin(pi/4)) 这样的数字精确到全精度,因为 MATLAB 通常会使用浮点运算计算调用 vpa 中的数字,因此只能精确到大约 16 位。

但是,这里似乎是正确的。例如,我们知道

sin(pi/4) == sqrt(2)/2

让我们测试一下这个结果。我将使用 100 位精度,比较 vpa 和我自己的 HPF 工具。

>> vpa(sin(pi/4),100)
ans =
0.7071067811865475244008443621048490392848359376884740365883398689953662392310535194251937671638207864

>> vpa(sqrt(sym(2))/2,100)
ans =
0.7071067811865475244008443621048490392848359376884740365883398689953662392310535194251937671638207864

>> sqrt(hpf(2,100))/2
ans =
0.7071067811865475244008443621048490392848359376884740365883398689953662392310535194251937671638207864

>> sin(hpf('pi',100)/4)
ans =
0.7071067811865475244008443621048490392848359376884740365883398689953662392310535194251937671638207864

所以,我的猜测是解析器已将输入识别为符号工具箱可以更准确地计算的内容。正如我之前所说,但要小心。 sin(pi/12)是什么?

>> vpa(sin(pi/12),100)
ans =
0.25881904510252073947640383266843855381011962890625

>> vpa('sin(pi/12)',100)
ans =
0.2588190451025207623488988376240483283490689013199305138140032073150569747488019969223679746942496655

>> vpa(sin(sym('pi')/12),100)
ans =
0.2588190451025207623488988376240483283490689013199305138140032073150569747488019969223679746942496655

>> sin(hpf('pi',100)/12)
ans =
0.2588190451025207623488988376240483283490689013199305138140032073150569747488019969223679746942496655

看到在第一种情况下,解析器没有拯救我们。在其他情况下,我强制 MATLAB 计算正确的值。事实上,稍加努力就能得到 sin(pi/12) 的值,如 sqrt(2)*(sqrt(3) - 1)/4。

>> DefaultNumberOfDigits 100
>> (sqrt(hpf(3)) - 1)*sqrt(hpf(2))/4
ans =
0.2588190451025207623488988376240483283490689013199305138140032073150569747488019969223679746942496655

关键是,不要相信解析器会把你救在这里。

编辑:作为对 Amro 评论的测试,我恭敬地声明 MATLAB 正在做一些有趣的事情。看到 vpa 能够返回正确的前 100 位 pi,即使将 pi 作为双精度数传递也是如此。由于 pi(作为双精度数)在第 16 位十进制数字之后是不正确的,因此发生了一些可疑的事情。

>> vpa(pi,100)
ans =
3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117068

>> vpa('pi',100)
ans =
3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117068

vpa('pi',100) - vpa(pi,100)
ans =
0.0

作为对这一事实的测试,让我们看看 HPF 发现了什么。 HPF 实际上采用 IEEE 754 值,以双精度形式存储,然后将其转换为 HPF 数。

>> hpf(pi,100)
ans =
3.141592653589793115997963468544185161590576171875

>> hpf('pi',100)
ans =
3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117068

>> hpf('pi',100) - hpf(pi,100)
ans =
0.0000000000000001224646799147353177226065932275001058209749445923078164062862089986280348253421170679821480800000000

很明显,MATLAB 能够将 pi 识别为不仅仅是将传入的双精度值。

编辑2:

事实上,一点玩意就告诉我这里发生了什么。 VPA 是一个棘手的问题,而不是解析器。考虑分数 7/13。如果我们将其构建为 double,然后打印出存储在其全部荣耀中的浮点值,我们会发现它并不完全准确。这符合预期。

>> sprintf('%.100f',7/13)
ans =
0.5384615384615384359179302009579259902238845825195312500000000000000000000000000000000000000000000000

7/13 是一个重复的十进制值。以下是正确的数字:

>> vpa('7/13',100)
ans =
0.5384615384615384615384615384615384615384615384615384615384615384615384615384615384615384615384615385

现在,假设我们尝试创建相同的数字。在这里,我会以双精度形式传递 7/13,但我会在底部的十进制数字中出错

>> sprintf('%.100f',0.538461538461538461777777777)
ans =
0.5384615384615384359179302009579259902238845825195312500000000000000000000000000000000000000000000000

在这里我们看到 vpa 捕获并纠正了我犯的“错误”,认识到我传入的值实际上与我在 7/13 传入时的值相同。

>> vpa(0.538461538461538461777777777,100)
ans =
0.5384615384615384615384615384615384615384615384615384615384615384615384615384615384615384615384615385

当然,如果我将值作为字符串传递,那么 vpa 会出错。

>> vpa('0.538461538461538461777777777',100)
ans =
0.538461538461538461777777777

这解释了为什么 vpa 能够捕获并正确计算 vpa(sin(pi/4),100),达到所要求的全部精度。 sin(pi/4) 被计算为双精度,但随后 vpa 将其视为与 sqrt(2)/2 的双精度版本相同的数字。

当然要小心。例如,vpa 不够聪明,无法捕捉到这个简单的 pi 变化。

>> vpa(pi + 1,100)
ans =
4.141592653589793115997963468544185161590576171875

>> vpa(pi,100)
ans =
3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117068

【讨论】:

  • 只是为了消除任何混淆,解析器不会将vpa(sin(pi/4)) 中的输入与调用任何其他函数myFcn(sin(pi/4)) 有任何不同... MATLAB 在传递结果之前首先以双精度计算表达式到功能。在有限精度下,将结果四舍五入到它可以表示的最接近的十进制值,VPA 将其视为精确值。到那时,它可以计算出尽可能多的小数点...我相信@AndrewJanke 在我链接到的问题中给出了很好的解释。
  • 非常感谢!非常感谢您的详细回答。
  • @woodchips:看来您也是正确的。当您使用数字表达式调用vpa 时,它在内部调用vpa(sym(expr))。现在,当使用数字调用 sym 时,它会尝试将数字转换为“有理”形式以补偿舍入错误(flag = 'r' 参数),这解释了为什么 VPA “纠正”了您所显示的错误.如果您想获得与 HPF 工具相同的结果,我想应该是:vpa(sym('pi')-sym(pi,'d'), 100)
  • @woodchips:难怪我错过了这个细节,查看previous versions 中的文档,这种转换并不明显,直到 R2012a 文档才提到这个事实很清楚...... +1 以获得真正有启发性的答案
【解决方案2】:

我不是 MatLab 专家,但如果没有引号,您会将 sqrt(5)/pi结果 传递给 vpa()

  vpa(sqrt(5)/pi)
= vpa(0.7117625434171772)

使用引号,您将 表达式 sqrt(5)/pi(未经计算且以精确形式)传递给 vpa(),然后告诉 MatLab 以可变精度计算 sqrt(5)/pi

【讨论】:

  • 谢谢!那讲得通。但是,为什么当我在vpa中简单地输入sin(pi/4)时,我是否使用引号并不重要?
  • 这应该很重要。正如@Ben Volgt 指出的那样,sin(pi / 4) = sin(45 deg) = sqrt(2) / 2,这是不合理的,不会用双精度浮点数精确逼近。
  • 读完the docs后听起来vpa返回的是双精度结果,这种情况下只有中间值的精度有所提高。
  • 太棒了。非常感激!我想我现在明白了:)
【解决方案3】:

如果你得到完全相同的答案,你就不需要变精度算术。

但是,sin(pi/4) 应该正好是sqrt(2)/2,这是不合理的。你不应该从不同的精度得到完全相同的答案。也许你应该检查你是如何显示(和四舍五入)结果的。

【讨论】:

  • 谢谢。我现在明白了:)。欣赏!
【解决方案4】:

numeric to symbolic 转换的最新文档有答案。

sym 尝试将浮点输入中的舍入误差纠正为 返回确切的符号形式。具体来说, sym 更正四舍五入 与 p/q、pπ/q、(p/q)^(1/2)、2^q、 和 10^q,其中 p 和 q 是中等大小的整数。

因此,sin(pi/4)2^(1/2)/2(1/2)^(1/2),因此 vpa 命令可以识别它。但是,根据文档,sqrt(5)/pi 不是可识别的输入表单。

【讨论】:

    猜你喜欢
    • 2014-06-24
    • 2013-09-06
    • 1970-01-01
    • 2010-11-16
    • 2011-02-10
    • 2012-08-23
    • 1970-01-01
    相关资源
    最近更新 更多