Tell me more ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.
up vote 10 down vote favorite
1

I am trying to plot a polynomial of degree 29 on the domain [0,1], with fairly large coefficients:

poly[z_] = -1.1126829840302355` + 113.28783661058498` z - 
9878.742379213338` z^2 + 584715.8646149524` z^3 - 
2.280647160113914`*^7 z^4 + 6.176933520283158`*^8 z^5 - 
1.217581378062843`*^10 z^6 + 1.811582263559531`*^11 z^7 - 
2.0920133095023196`*^12 z^8 + 1.9158620445090305`*^13 z^9 - 
1.414955836879161`*^14 z^10 + 8.539129365981781`*^14 z^11 - 
4.254563022912091`*^15 z^12 + 1.764182366816184`*^16 z^13 - 
6.125080435776876`*^16 z^14 + 1.7883504482275766`*^17 z^15 - 
4.403320010637656`*^17 z^16 + 9.154951756734264`*^17 z^17 - 
1.6068672087698447`*^18 z^18 + 2.3765393965161196`*^18 z^19 - 
2.950846281328122`*^18 z^20 + 3.0579497598096415`*^18 z^21 - 
2.6220913470110597`*^18 z^22 + 1.837651151556163`*^18 z^23 - 
1.0344252684666292`*^18 z^24 + 4.5602474296077024`*^17 z^25 - 
1.5155510563521117`*^17 z^26 + 3.568596763872067`*^16 z^27 - 
5.304183668348243`*^15 z^28 + 3.7404713997980006`*^14 z^29

The problem is when I try to plot the polynomial all seems fine near 0 but then mid way through the the unit interval I observe some erratic behaviour as shown below. Can anyone tell me why this is and how I can avoid it?

Plot[{poly[z]}, {z, 0, 1}, WorkingPrecision -> precision]

I would expect the polynomial to be a smooth curve

share|improve this question
5  
Have you considered that you're working with a 29th degree polynomial with insanely huge coefficients? It's not surprising to me that you get a plot like that. It may seem "nice" early on, but that's because z^29 (and other high-order terms) is very small near 0. Near approximately 0.85, it rapidly approaches 1. – Mike Bantegui Mar 18 '12 at 2:14
3  
@Mithc, I will have to disagree with you. Just because the terms of the ploynomial get large as z goes to 1 does not mean you should expect this type of erratic behaviour. Indeed the polynomial can have no more than 29 real roots! So why you expect a plot like this is beyond me. You should not have voted the question down. I suspect a numerical issue here and was hoping somebody could help me pin it down. – aukie Mar 18 '12 at 3:09
6  
It's cancellation error. Per advice in the response, use higher precision and/or rationalize the coefficients. And be sure you know what you are doing in terms of what might comprise an "expected" outcome, because the result might not mean much if the actual coefficients really have error intervals associated with them. – Daniel Lichtblau Mar 19 '12 at 1:16

migrated from stackoverflow.com Mar 18 '12 at 20:11

1 Answer

up vote 18 down vote accepted

If you Rationalize your real numbers you will be able to use Mathematica's arbitrary precision engine:

poly2 = Rationalize[poly[z], 0];

Plot[Evaluate[poly2], {z, 0, 1}, WorkingPrecision -> 50]

Mathematica graphics

Mathematica has two kinds of numeric calculations: machine precision, and arbitrary precision. Machine precision is fast but is limited to ~16 digits (53 binary digits) and may lose precision during a calculation. Mathematica also does not track the precision of a result.

Numbers entered as 1.234 or 1234` are taken to be machine precision.

Arbitrary precision is slower, but Mathematica will track the precision of calculations, often using more calculations as necessary to preserve precision, and print the precision of the result.

Exact values such as 1234 and 1/2 can be used in arbitrary precision calculations. Numbers can also be entered with e.g. 1.234`20 specifying 20 digits of precision, and these will automatically use arbitrary precision if all other values are either exact or arbitrary.

Precision can be checked with the function Precision:

Precision /@ {1.234, 1234`, 1.234`20, 7}

{MachinePrecision, MachinePrecision, 20., ∞ }

Precision can only be preserved if all values in a calculation have at least that precision. Also, arbitrary precision arithmetic may be used with numbers having a precision less than MachinePrecision -- Mathematica will show the true precision of the result.

Precision[1.234 + 7]

Precision[1.234`20 + 1.234`12]

MachinePrecision

12.301

Precision can be set with SetPrecision. It is probably better to use this rather than Rationalize to put numbers into a form that the arbitrary precision engine will use, because the latter will be manufacturing false precision.

Applying this to your problem:

poly3 = SetPrecision[poly[z], 15];

Plot[Evaluate[poly3], {z, 0, 1}, WorkingPrecision -> 50]

During evaluation of In[91]:= Plot::precw: The precision of the argument function <<>> is less than WorkingPrecision (50.`). >>

This is an important warning because it lets you know that your results may not be valid.

See this tutorial for more information about precision. Take time to understand the difference between Mathematica's meanings of Accuracy and Precision.

share|improve this answer
@Wizard, Thank you. If possible can you give me some insight into the problem and the solution. – aukie Mar 18 '12 at 3:22
@aukie see the updated answer – Mr.Wizard♦ Mar 18 '12 at 12:05
4  
To complement this, calculating the polynomial value in Horner form (HornerForm) increases precision somewhat, but not enough. – Szabolcs Mar 18 '12 at 22:26

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.