Re: [R] Solving 100th order equation

From: Hans W. Borchers <hwborchers_at_gmail.com>
Date: Sun, 25 May 2008 22:38:25 +0000 (UTC)

Rolf Turner <r.turner <at> auckland.ac.nz> writes:

>
> [...]
>
> However uniroot(p,c(1.995,2.005)) gives
>
> $root
> [1] 1.999993
>
> $f.root

> [1] -4.570875e+24
>
> $iter
> [1] 4
>
> $estim.prec
> [1] 6.103516e-05
>
> What a difference 7.214144e-06 makes! When you're dealing with
> polynomials of degree 100.
>

I don't think R is the right tool to solve this kind of questions that belong to the realm of Computer Algebra Systems. Yacas is much to weak for such high-order polybomials, but we can apply more powerful CASystems, for instance the free Maxima system. Applying

(%i) nroots(x^100 - 2*x^99 + 10*x^50 + 6*x - 4000, minf, inf)

immediately shows that there are only *two* real solutions, and then

(%i) a: realroots(x^100 - 2*x^99 + 10*x^50 + 6*x - 4000, 1e-15); (%i) float(a)

(%o) [x=-1.074126672042147,x=1.999999999999982]

will provide real solutions with 15 decimals (does not change when more decimals are used). So the difference that counts is actually much smaller.

The complex roots -- if needed -- will require special treatment.

I believe it would be fair to point such questions to Computer Algebra mailing lists and not try to give the appearance they could be solved satisfyingly with R. The same as a complicated statistics question in a Matlab or Mathematica mailing list should be pointed to R-help.

Hans Werner

> Bozhemoi!
>
> cheers,
>
> Rolf Turner



R-help_at_r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Sun 25 May 2008 - 22:42:10 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Mon 26 May 2008 - 00:30:40 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.

list of date sections of archive