Main Page | Report this Page
 
   
Science Forum Index  »  Math - Symbolic Forum  »  Maxima - precision
Page 1 of 1    
Author Message
adam majewski
Posted: Sun May 20, 2007 2:59 am
Guest
Hi

I'm computing complex roots of polynomial:
P(n):=if n=0 then 0 else P(n-1)^2+c;
bigfloat(allroots(%i*P(3)));

How can I set precision ( number of significant digits) ?
Changing fpprec seems to have no effect .

Regards

Adam


===========================================================

(%i1) P(n):=if n=0 then 0 else P(n-1)^2+c;
(%o1) P(n):=if n=0 then 0 else P(n-1)^2+c
(%i2) bigfloat(allroots(%i*P(3)));
(%o2)
bigfloat([c=0.0,c=0.74486176661974*%i-0.12256116687665,c=-0.74486176661974*%i-0.12256116687665,c=-
1.754877666246693])
(%i3) fpprec : 100;
(%o3) 100
(%i4) bigfloat(allroots(%i*P(3)));
(%o4)
bigfloat([c=0.0,c=0.74486176661974*%i-0.12256116687665,c=-0.74486176661974*%i-0.12256116687665,c=-
1.754877666246693])
(%i5) allroots(%i*P(3));
(%o5)
[c=0.0,c=0.74486176661974*%i-0.12256116687665,c=-0.74486176661974*%i-0.12256116687665,c=-
1.754877666246693]
(%i6) allroots(bigfloat(%i*P(3)));
`allroots': not a polynomial:
bigfloat(%i*((c^2+c)^2+c))
-- an error. To debug this try debugmode(true);
(%i7) fpprec : 2;
(%o7) 2
(%iCool allroots(bigfloat(%i*P(3)));
`allroots': not a polynomial:
bigfloat(%i*((c^2+c)^2+c))
-- an error. To debug this try debugmode(true);
(%i9) bigfloat(allroots(%i*P(3)));
(%o9)
bigfloat([c=0.0,c=0.74486176661974*%i-0.12256116687665,c=-0.74486176661974*%i-0.12256116687665,c=-
1.754877666246693])
(%i10) realroots(%);
(%o10)
[bigfloat([c=0.0,c=0.74486176661974*%i-0.12256116687665,c=-0.74486176661974*%i-0.12256116687665,c=
-1.754877666246693])=0]
(%i11) fpprec : 2;
(%o11) 2
(%i12) realroots(P(3));
(%o12) [c=-58883923/33554432,c=0]
(%i13) float(%), numer;
(%o13) [c=-1.754877656698227,c=0.0]
(%i14)
rjf
Posted: Sun May 20, 2007 7:11 pm
Guest
On May 20, 12:59 am, adam majewski <adamm...@o2.pl> wrote:
Quote:
Hi

I'm computing complex roots of polynomial:
P(n):=if n=0 then 0 else P(n-1)^2+c;
bigfloat(allroots(%i*P(3)));

How can I set precision ( number of significant digits) ?
Changing fpprec seems to have no effect .

You are computing the roots of a polynomial in machine precision.
Allroots is just a transcription of the Jenkins-Traub floating point
polynomial zero-finder, and does not use bigfloats at all.

given the machine-float roots, you then convert them to bigfloat
precision. This does not make them any more accurate. At least you
would do this if you were to spell bfloat(...) correctly.

What you want is a program different from allroots that computes the
roots to a precision that depends on fpprec. This is not, so far as I
know, part of the Maxima system, though a version could easily be
produced.
It would probably have some glitches because zero-finding programs
have to deal with nasty problems.
But maybe your polynomials are not nasty.

What you can do is choose a zero computed from allroots, and then
refine it.
Usually you can use newton iteration, and run that until you are happy
with the answer.
adam majewski
Posted: Sat Jun 02, 2007 11:18 am
Guest
rjf napisał(a):
Quote:
On May 20, 12:59 am, adam majewski <adamm...@o2.pl> wrote:
Hi

I'm computing complex roots of polynomial:
P(n):=if n=0 then 0 else P(n-1)^2+c;
bigfloat(allroots(%i*P(3)));

How can I set precision ( number of significant digits) ?
Changing fpprec seems to have no effect .

You are computing the roots of a polynomial in machine precision.
Allroots is just a transcription of the Jenkins-Traub floating point
polynomial zero-finder, and does not use bigfloats at all.

given the machine-float roots, you then convert them to bigfloat
precision. This does not make them any more accurate. At least you
would do this if you were to spell bfloat(...) correctly.

What you want is a program different from allroots that computes the
roots to a precision that depends on fpprec. This is not, so far as I
know, part of the Maxima system, though a version could easily be
produced.
It would probably have some glitches because zero-finding programs
have to deal with nasty problems.
But maybe your polynomials are not nasty.

What you can do is choose a zero computed from allroots, and then
refine it.
Usually you can use newton iteration, and run that until you are happy
with the answer.


Thx for answer. If I properly understand :


P(n):=if n=0 then 0 else P(n-1)^2+c;
load(newton1);
clist:allroots(%i*P(6));
a:rhs(clist[4]);
(%i9) newton(%i*P(3),c,a,0.000000000000000000001);Maxima encountered a
Lisp error: Error in PROGN [or a callee]: The storage for CONS is
exhausted.Currently, 82780 pages are allocated.Use ALLOCATE to expand
the space.Automatically continuing.To reenable the Lisp debugger set
*debugger-hook* to nil.


Adam
rjf
Posted: Sat Jun 02, 2007 4:59 pm
Guest
Newton iteration does not always converge. Why should yours converge?
RJF

adam majewski wrote:
Quote:
rjf napisał(a):
On May 20, 12:59 am, adam majewski <adamm...@o2.pl> wrote:
Hi

I'm computing complex roots of polynomial:
P(n):=if n=0 then 0 else P(n-1)^2+c;
bigfloat(allroots(%i*P(3)));

How can I set precision ( number of significant digits) ?
Changing fpprec seems to have no effect .

You are computing the roots of a polynomial in machine precision.
Allroots is just a transcription of the Jenkins-Traub floating point
polynomial zero-finder, and does not use bigfloats at all.

given the machine-float roots, you then convert them to bigfloat
precision. This does not make them any more accurate. At least you
would do this if you were to spell bfloat(...) correctly.

What you want is a program different from allroots that computes the
roots to a precision that depends on fpprec. This is not, so far as I
know, part of the Maxima system, though a version could easily be
produced.
It would probably have some glitches because zero-finding programs
have to deal with nasty problems.
But maybe your polynomials are not nasty.

What you can do is choose a zero computed from allroots, and then
refine it.
Usually you can use newton iteration, and run that until you are happy
with the answer.


Thx for answer. If I properly understand :

P(n):=if n=0 then 0 else P(n-1)^2+c;
load(newton1);
clist:allroots(%i*P(6));
a:rhs(clist[4]);
(%i9) newton(%i*P(3),c,a,0.000000000000000000001);Maxima encountered a
Lisp error: Error in PROGN [or a callee]: The storage for CONS is
exhausted.Currently, 82780 pages are allocated.Use ALLOCATE to expand
the space.Automatically continuing.To reenable the Lisp debugger set
*debugger-hook* to nil.


Adam
 
Page 1 of 1       All times are GMT - 5 Hours
The time now is Thu Jan 08, 2009 12:00 am