这是mathematica中的NSolve中的错误吗？

Question

人们会期望并希望如果你要求Mathematica找到多项式的根,它应该给出相同的(近似的)答案,无论你是否象征性地做这个,然后找到这些确切答案的数值近似值,或者你是否用数字表示.这是一个例子(在Mathematica 7OS X上运行),这种情况严重失败:

poly = -112 + 1/q^28 + 1/q^26 - 1/q^24 - 6/q^22 - 14/q^20 - 25/q^18 -
  38/q^16 - 52/q^14 - 67/q^12 - 81/q^10 - 93/q^8 - 102/q^6 - 108/
  q^4 - 111/q^2 - 111 q^2 - 108 q^4 - 102 q^6 - 93 q^8 - 81 q^10 -
  67 q^12 - 52 q^14 - 38 q^16 - 25 q^18 - 14 q^20 - 6 q^22 - q^24 +
  q^26 + q^28;

Total[q^4 /. NSolve[poly == 0, q]] - Total[q^4 /. N[Solve[poly == 0, q]]]

(注意:这实际上是一个Laurent多项式,如果你乘以大q问题的力量就会消失.)

这里的最后一行只是一个证明,所发现的解决方案非常不同; 事实上,这是我们在努力解决的问题中尝试计算的数量.

如果你的输出仔细观察NSolve[poly == 0, q]和N[Solve[poly == 0, q],你会看到NSolve只给出了54根,而不是预期56.并不只是错过了重复的根或任何东西; 它错过了两个最大的根(大约+/- 1.59)

这是Mathematica中的错误吗？有没有人解释为什么会这样？

Answer 1

不幸的是,你的期望没有根据.

如你所说,Solve[]给出一个精确的解决方案并N[]引入一个小错误,但只有一次,当你评估它时.NSolve[]另一方面,使用来自get get的数值近似,因此可以累积显着的舍入误差.

您还受到计算的默认精度的限制,这可能导致数值方法完全失败,例如缺少根(参见Wilkinson多项式).您可以通过告诉NSolve[]使用更高的精度来抵消这种情况:

In[1]  := Total[q^4 /. NSolve[poly == 0, q, WorkingPrecision -> 50]] - 
          Total[q^4 /. N[Solve[poly == 0, q]]]
Out[1] := 0. - 3.66374*10^-15 I

使用数值方法时,记住错误始终很重要.因为对于从解决长多项式到大矩阵对角化到奇怪函数等的集成等各种各样的数值分析问题都是如此.没有一种正确的方法,并且Mathematica需要被告知,例如,提高WorkingPrecision ,或应用不同的数字技术.