On Thu, 7 May 2009, Oscar wrote:
< Hi all,
< I want to get the roots of a polynomial. I used the allroot funcion
< (as I will show with a code below), and then I mapped the same
< polynomial using the value of a root returned from allroots, but it
< didn't work.
<
< (%i200) g(x) := x^5-x^3-10^(-15);
< allroots(%i*g(x));
< (%o200) g(x):=x^5-x^3-10^(-15)
< (%o201) [x=8.6602540378443868*10^-6*%i+4.9999999996666667*10^-6,x=-4.5460856822775728*10^-22*%i-1.0000000000333334*10^-5,x=4.9999999996666658*10^-6-8.6602540378443851*10^-6*%i,x=1.0-1.8541339457159821*10^-17*%i,x=1.8539432108527007*10^-17*%i-1.0]
<
< (%i202) map(g, [8.6602540378443868*10^-6*%i+4.9999999996666667*10^-6]);
< (%o202) [(8.6602540378443868*10^-6*%i+4.9999999996666667*10^-6)^5-(8.6602540378443868*10^-6*%i+4.9999999996666667*10^-6)^3-1/1000000000000000]
You should apply `expand' to this. It should simplify to something quite
small. It is also a good idea to `reset();' a session if you are having
problems -- start calculations with a clean piece of paper.
(%i3) g(x) := x^5-x^3-10^(-15)$
(%i4) r : allroots(g(x));
(%o4) [x = -1.0000000000333335E-5,
x = 8.660254037844386E-6*%i+4.999999999666667E-6,
x = 4.999999999666667E-6-8.660254037844386E-6*%i,
x = -.9999999999999994,x = 1.0]
(%i5) map(lambda([x],expand(g(rhs(x)))), r);
(%o5)
[5.916456789157588E-31,2.368316254054144E-32-3.055051027801933E-31*%i,
3.055051027801933E-31*%i+2.368316254054144E-32,1.102230246251565E-16,
-1.118215802998748E-16]
(%i6) map(lambda([x],expand(abs(g(rhs(x))))), r);
(%o6)
[5.916456789157588E-31,3.064217029073791E-31,3.064217029073791E-31,
1.102230246251565E-16,1.118215802998748E-16]
You can see that the modulus of g(r) is < 1e-15 in all cases.
Here is the example, but with an errant %i in front. Perhaps the
developers can comment on why allroots finds slightly different roots. I
would have thought that it would have first normalised %i*g(x) to a
monic polynomial (i.e. g(x)) before root finding.
(%i7) r : allroots(%i*g(x));
(%o7) [x = 8.660254037844386E-6*%i+4.999999999666666E-6,
x = -4.546085682277573E-22*%i-1.0000000000333333E-5,
x = 4.9999999996666655E-6-8.660254037844385E-6*%i,
x = 1.0-1.854133945715982E-17*%i,
x = 1.853943210852701E-17*%i-.9999999999999994]
(%i8) map(lambda([x],expand(abs(g(rhs(x))))), r);
(%o8) [2.04548464889026E-31,2.026200846528237E-31,5.680450389142161E-31,
2.25119816213002E-16,1.170503900707085E-16]
Leo
<
< The map of a root has to be zero, doesn't it?
< I have already tried with allroots(g(X)), but I got some strange
< results (can someone explain me why?)
<
< By the way: my work with maxima is to find some polynomials and value
< to which some rounding error happends, can someone give me some source
< in the literature the approachs this topic (with examples)? I want to
< identify interesting polynomials to test my new benchmark for floating
< point calculation, which evaluates the quality of the result
< (accuracy). I hope someone can help me.
<
< Thank you,
< Oscar
< _______________________________________________
< Maxima mailing list
< Maxima at math.utexas.edu
< http://www.math.utexas.edu/mailman/listinfo/maxima
<
<
--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.