1
$\begingroup$

Background: I am trying to compute the Lyapunov exponent for the logistic map given by \begin{equation*} x_{n+1}=F(x_n)=ax_n(1-x_n) \end{equation*} The formula for the Lyapunov exponent for a given initial condition $x_0$ is given as \begin{equation*} \lambda(x_0) = \lim_{n\rightarrow\infty}\frac{1}{n}\ln\left|{\frac{dF^n(x_0)}{dx_0}}\right| \end{equation*} While I do not intend to check for large values of $n$, as the degree of $F^n$ increases as $2^n$, I would like to check up to $n=10$ or so. I would also like to plot $\lambda$ vs $n$ (removing the limit) but that part is not included in the code.

My entire code is as follows ( a, nmax are taken as inputs earlier):

forwardList = NestList[ a # (1 - #)&, x, nmax];
derivatives = Table[ D[ forwardList[[2]], x], {i, nmax}]
Do[derivatives[[i]] = a (1 - 2forwardList[[i]]) derivatives[[i - 1]], {i, 2, nmax}]
Log[ Abs[ derivatives]] /. x -> 0.3

I did not divide by $n$ just for convenience.

Problem: Once I tried using Collect for checking the numbers and I got different answers.

Consider the following code snippet:

forwardList = Collect[NestList[ a # (1 - #)&, x, nmax], x];
Collect[forwardList - NestList[ a # (1 - #)&, x, nmax], x]

I noted that the latter collect is giving me a non-zero list for list elements after the fourth one or so. Since the logistic map is chaotic, making precise calculations is very important. An example: for a = 3.9, nmax = 10 the last element of forwardList turns out to be ~ 5.5.
Instead of we replace the definition of forwardList (the topmost one) as follows, we get ~ 809.

forwardList = Collect[ NestList[ a # (1 - #)&, x, nmax], x];

Which of the answers is correct and why? Does Collect work properly on lists?

A correct alternative implementation would be welcome. My current algorithm uses a lot of RAM :( so I think there must be a better way.

$\endgroup$

1 Answer 1

2
$\begingroup$

It seems to be because of numerical errors. There are A LOT of operations so there is plenty of opportunity for errors to grow.

To get around this you can use exact arithmetic which Mathematica does if it is given exact numbers, for instance a=3.9 is stored as a floating point number which is inexact causing all computation done to be inexact too. a=39/10 is exact and all computation will be done with exact arithmetic.

a=39/10;
nmax=10;
forwardList = Collect[NestList[a*#*(1 - #) &, x, nmax], x];
Collect[forwardList - NestList[a*#*(1 - #) &, x, nmax], x]
(* {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0} *)

In these cases you usually want to put in numerical values as late as possible so you can do the main computation with a unassigned once for each n you are interested in and then replace a with any value you want.

When you have inexact a or x you can use Rationalize to get something close that is exact. And most functions have a WorkingPrecision option that can help out.

lyapunov

Click Learning Resources in the Numerical Evaluation And Precision Guide for an overview of how all this works.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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