Re: Bug in NSolve?
- To: mathgroup at smc.vnet.net
- Subject: [mg33843] Re: [mg33788] Bug in NSolve?
- From: "Fred Simons" <f.h.simons at tue.nl>
- Date: Fri, 19 Apr 2002 02:29:46 -0400 (EDT)
- References: <200204160750.DAA20077@smc.vnet.net>
- Sender: owner-wri-mathgroup at wolfram.com
It looks like a bug, maybe related to some other strange behaviour I will
show in the examples. Maybe the problem has to do with small but not very
small imaginary parts of the solutions
We define a polynomial of degree depneding on two parameters such that the
imaginary parts of the two roots are-1/b and 1/b.
In[1]:=
Clear[a, b];
pol[a_, b_] = 1 + a^2 - 2*a*b*x + b^2*x^2;
Solve[pol[a, b] == 0, x]
Out[3]=
{{x -> (-I + a)/b}, {x -> (I + a)/b}}
In[4]:=
a = 45654; b = 1;
NSolve[pol[a, b] == 0, x]
pol[a, b] /. %
Out[5]=
{{x -> 45654. - 0.999999957203608*I},
{x -> 45654. + 0.999999957203608*I}}
Out[6]=
{0. + 0.*I, 0. + 0.*I}
This is OK. But with one more digit in a we get incorrect results:
In[7]:=
a = 456544; b = 1;
NSolve[pol[a, b] == 0, x]
pol[a, b] /. %
Out[8]=
{{x -> 456543.9999666127}, {x -> 456544.0000333873}}
Out[9]=
{1.000030517578125, 0.999969482421875}
Fortunately, using arbitrary precision numbers helps:
In[10]:=
a = 456544; b = 1;
NSolve[pol[a, b] == 0, x, 20]
pol[a, b] /. %
Out[11]=
{{x -> 456543.999999999999999999771731`14.0141 -
0.9999999999999999958`8.3546*I},
{x -> 456543.999999999999999999771731`14.0141 +
0.9999999999999999958`8.3546*I}}
Out[12]=
{-1.14844`0*^-9 + 9.2294`0*^-10*I,
-1.14844`0*^-9 + -9.2294`0*^-10*I}
Now we look at b:
In[13]:=
a = 456544; b = 276355;
NSolve[pol[a, b] == 0, x, 20]
pol[a, b] /. %
Out[14]=
{{x -> 1.6520200466790903005192583234`14.014 -
3.618534131823198421`8.3545*^-6*I},
{x -> 1.6520200466790903005192583234`14.014 +
3.618534131823198421`8.3545*^-6*I}}
Out[15]=
{-1.271611`0*^-8 + 9.90684`0*^-9*I,
-1.271611`0*^-8 + -9.90684`0*^-9*I}
One digit more in b gives the strange and incorrect result observed by
rodgerro at library2.airnews.net :
In[16]:=
a = 456544; b = 2763556;
NSolve[pol[a, b] == 0, x, 20]
pol[a, b] /. %
Out[17]=
{{x -> 0.1652016459952322297793120673`14.0137}}
Out[18]=
{1.`2.0927}
One more increasing precision helps:
In[19]:=
a = 456544; b = 2763556;
N[NSolve[pol[a, b] == 0, x, 39]]
pol[a, b] /. %
Out[20]=
{{x -> 0.16520164599523224 - 3.6185262755666976*^-7*
I}, {x -> 0.16520164599523224 +
3.6185262755666976*^-7*I}}
Out[21]=
{0. + 0.*I, 0. + 0.*I}
One digit more in b yields incorrect results:
In[22]:=
a = 456544; b = 27635567;
N[NSolve[pol[a, b] == 0, x, 39]]
pol[a, b] /. %
Out[23]=
{{x -> 0.016520160415018807}}
Out[24]=
{0.99993896484375}
And so on?
Fred Simons
Eindhoven University of Technology
- References:
- Bug in NSolve?
- From: rodgerroSPAMNOT.siteconnect.com@library2.airnews.net
- Bug in NSolve?