0

我在matlab中实现一个代码来求解二次方程，使用解析公式：

在此处输入图像描述

这是代码：

clear all 
format short 
a=1; b=30000000.001; c=1/4; 
rdelta=sqrt(b^2-4*a*c); 
x1=(-b+rdelta)/(2*a);
x2=(-b-rdelta)/(2*a);
fprintf(' Roots of the polynomial %5.3f x^2 + %5.3f x+%5.3f \n',a,b,c)  
fprintf ('x1= %e\n',x1)
fprintf ('x2= %e\n\n',x2)
valor_real_x1= -8.3333e-009;
valor_real_x2= -2.6844e+007;

error_abs_x1 = abs (valor_real_x1-x1);
error_abs_x2 = abs (valor_real_x2-x2);

error_rel_x1 = abs (error_abs_x1/valor_real_x1);
error_rel_x2 = abs (error_abs_x2/valor_real_x2);

fprintf(' absolute_errorx1 = |real value - obtained value| = |%e - %e| = %e \n',valor_real_x1,x1,error_abs_x1)  
fprintf(' absolute_errorx2 = |real value - obtained value| = |%e - %e| = %e \n\n',valor_real_x2,x2,error_abs_x2) 

fprintf(' relative error_x1 = |absolut error / real value| = |%e / %e| = %e \n',error_abs_x1,valor_real_x1,error_rel_x1 ) 
 fprintf(' relative_error_x2 = |absolut error / real value|  = |%e / %e| = %e \n',error_abs_x2,valor_real_x2,error_rel_x2)

我遇到的问题是它给了我一个精确的解决方案，即对于值a = 1，b = 30000000,001 c = 1/4，根的值是：

Roots of the polynomial 1.000 x^2 + 30000000.001 x+0.250 
 x1= -9.313226e-009
 x2= -3.000000e+007

知道多项式根的精确值为：

x1= -8.3333e-009
x2= -2.6844e+007

这给了我在计算的绝对和相对精度方面的以下错误：

 absolute_errorx1 = |real value - obtained value| = |-8.333300e-009 - -9.313226e-009| = 9.799257e-010 
 absolute_errorx2 = |real value - obtained value| = |-2.684400e+007 - -3.000000e+007| = 3.156000e+006 

 relative error_x1 = |absolut error / real value| = |9.799257e-010 / -8.333300e-009| = 1.175916e-001 
 relative_error_x2 = |absolut error / real value|  = |3.156000e+006 / -2.684400e+007| = 1.175682e-001

我的问题是：是否有获得二次方程根的最佳方法？即我可以更改我的代码以减少预期解决方案和结果解决方案之间的相对误差？

4

4 回答 4

7

在这种情况下直接使用二次公式会因减去两个幅度非常相似的值而导致数值精度的大量损失。这是因为表达式

sqrt(b*b - 4*a*c)

与 b 几乎相同。所以你应该只使用这两个根中的一个，一个不涉及减去两个非常接近的值，而对于另一个根，你可以使用（例如）二次根的乘积是 c/a 的事实。我会让你填补空白。

于 2012-04-04T13:53:26.947 回答
3

为什么这听起来像是数值分析第一堂课的家庭作业？

自从我年轻时已经有一段时间了，但我记得有一个技巧。无论如何，你错了。该多项式的真根是

solve('x^2 + 30000000.001*x + 0.25')
ans =
          -30000000.000999991666666666944442
 -0.0000000083333333330555578703796293981491

根在这里做得如何？

p = [1 30000000.001 1/4];
format long g
roots(p)
ans =
             -30000000.001
     -8.33333333305556e-09

这实际上看起来还不错。HPF 是怎么做的？

DefaultNumberOfDigits 64
a = hpf(1);
b = hpf('30000000.001');
c = hpf('0.25');

r1 = (-b + sqrt(b*b - 4*a*c))/(2*a)
r1 =
-0.000000008333333333055557870379629398149125529835186899898569329967

r2 = (-b - sqrt(b*b - 4*a*c))/(2*a)
r2 =
-30000000.000999991666666666944442129620370601850874470164813100101

是的，HPF 也可以很好地工作。

那么当你使用双精度数字和标准公式时会发生什么？是的，克拉波拉来了。

a = 1;
b = 30000000.001;
c = 0.25;

>> r1 = (-b + sqrt(b*b - 4*a*c))/(2*a)
r1 =
     -7.45058059692383e-09

>> r2 = (-b - sqrt(b*b - 4*a*c))/(2*a)
r2 =
             -30000000.001

同样，大量的减法抵消会侵蚀结果。（我似乎记得那是您在上一个问题中遇到的问题。）

您可以使用一个技巧。看到大解被很好地估计了，而不是接近于零的解。那么，如果您使用二次公式求解 Fliplr(p) 的根，会发生什么情况？这如何解决您的问题？当你这样做时，隐含地做了什么转换？（对不起，我不会做你的功课。我认为上面的提示已经足够了。）

于 2012-04-04T13:46:33.573 回答
1

我认为您的“真实”值可能是错误的（或者可能是精确的事情......我不知道）

a*(valor_real_x1^2)+b*(valor_real_x1)+c

ans =

   9.9999e-07

a*(valor_real_x2^2)+b*(valor_real_x2)+c

ans =

  -8.4720e+13
于 2012-04-04T13:41:03.260 回答
1

这个问题的一个很好的公式：

var q = sqrt(c*a)/b;
var f = .5 + .5 *sqrt(1-4*q*q);
var x1=-b*f/a;
var x2=-c/(f*b);
于 2015-11-03T15:44:03.690 回答