I am trying to make some simulations of chaotic systems, for instance X(k) = 4 X(k) (1 - X(k-1)) but I noticed that for all these systems, the loss of precision propagates exponentially, to the point that after 50 iterations, all values generated are completely wrong. I wrote some code in Perl using the BigNum library (providing hundreds of decimals accuracy) and it shows how dramatic standard arithmetic fails in this context.

You can check out the context, my code, and an Excel spreadsheet that illustrates the issue, here.

I am looking for a piece of code in Python that could nicely do the job, probably using some kind of BigNum library in Python? Anyone can make recommendations, or re-write my code in Python? Alternatively, how could this be done in R?

Thank you!

Source for picture: click here.

For arbitrary precision in many programming languages, check out this reference. Not sure if it is up-to-date and correct, but could not find anything about R. Some of these packages are not truly "arbitrary precision." More on this (for Python) here.