Re: Possible bug in NSolve[equation, variable, precission]
- To: mathgroup at smc.vnet.net
- Subject: [mg74394] Re: [mg74299] Possible bug in NSolve[equation, variable, precission]
- From: Daniel Lichtblau <danl at wolfram.com>
- Date: Mon, 19 Mar 2007 22:06:23 -0500 (EST)
- References: <200703170710.CAA19095@smc.vnet.net>
Julian Aguirre wrote: > Dear group, > > Mathematica 5.2 chokes solving numerically a polynomial equation. > > In[1] := $Version > Out[1]= 5.2 for Mac OS X (64 bit) (June 20, 2005) > > In[2]:= poly=171142046150220198693105489-16023210221608713837587916 > x-2020825892011586434364754 x^2+190894692033395024364972 > x^3+6039743423966949379761 x^4-568929229651998950400 > x^5-470066550477520896 x^6+2821109907456 x^7; > > In[3]:= poly2=Expand[poly/9]; > > In[4]:= NSolve[poly==0,x] > Out[4]= {-1211.83, -13.0015, -13.0014, 11.923, 12.0809, 12.2509, > 167826.} > > (* Up to this moment, everything is O.K. But *) > > In[5]:= NSolve[poly==0,x,20] > Out[5]= $Aborted (* after a loooong time *) > > (* However, the following works as expected*) > > In[6]:= x/.NSolve[poly2==0,x,20] > Out[6]= {-1211.8267955098487289, -13.001455891126, -13.001441554521, > 11.92303189062617, 12.08089051352363, 12.25087466630727, > 167826.26017849924816} > > Let me say that I have used Mathematica to solve thousands (probably > millions) of equations like the one above. There must be some magic in > the coefficients! > > Julian Aguirre > University of the Basque Country > Yes, it has to do with the coefficients and some other features of the polynomial. Adam Strzebonski and I investigated a similar set of examples reported to us around two weeks ago. The user who brought it to our attention also apparently does massive numbers of these in the course of a large simulation, estimateds to be in the tens or hundreds of thousands over a period of several hours. Over a period of many such runs a few such examples were culled. So it does indeed seem to be a bit of a rare phenomenon, though unfortunately not impossible to encounter. What we know at present is this appears to involve a long-standing problem in the numeric root-finding code used in the Mathematica kernel. The trouble seems to correlate with near-multiplicity (but not with actual multiplicity to within the requested bignum precision), along with bad scaling. Your example also indicates it may have something to do with the absolute magnitude of the coefficients, since rescaling by a factor of 9 (the gcd of the coefficients) makes the problem disappear in that case. At this time we do not know in more detail what is behind the hang. One thing to realize is that NSolve has dedicated special case handling for a single univariate polynomial. Thus one workaround is to add an equation in another variable e.g. y==x, to force NSolve to use a different approach. Other possibilities include rescaling, or setting up a companion matrix and finding the eigenvalues at whatever precision is of interest (this uses Lapack methods rather than a Jenkins-Traub based rootfinder). Perhaps the main issue, for practical purposes, is how you can tell when you will need a workaround. The best we have to offer is to use TimeConstrained on your NSolve, and in cases where it aborts, go to plan b. Daniel Lichtblau Wolfram Research
- References:
- Possible bug in NSolve[equation, variable, precission]
- From: "Julian Aguirre" <julian.aguirre@ehu.es>
- Possible bug in NSolve[equation, variable, precission]