# Calculating the Feigenbaum Constants

I would like to calculate the $δ$ and the $α$ Feigenbaum constants for the Logistic Map:

$$x_{n+1}=4μx_n(1-x_n)=f(x_n)$$ as part of an undergraduate assigment of mine. I have managed so far to calculate explicitly the 2-cycle and the value of $\mu$ for which it becomes unstable, leading to the bifurcation of the 4-cycle. My attempt of calculating the $δ$, which is defined as $$δ=\lim_{n \to \infty}\frac{μ_{n}-μ_{n+1}}{μ_{n+1}-μ_{n+2}}$$ is given below

f[x_Real, mu_Real, n_Real] := Block[{y = x}, Do[y = 4 mu y (1 - y), {n}]; y]
mu /. FindRoot[f[1/2, mu, 2^2] == 1/2, {mu, 0.9}]
delta[n_] := (mu[n - 1] - mu[n - 2])/(mu[n] - mu[n - 1])
mu[n_] := mu[n] == mu /. FindRoot[f[1/2, mu, 2^n] ==
1/2, {mu, {mu[n - 1] + (mu[n - 1] - mu[n - 2])/delta[n - 1],
mu[n - 1] + (mu[n - 1] - mu[n - 2])/(delta[n - 1] + .01)}},
AccuracyGoal -> 10]
Table[mu[n], {n, 0, 7}]
Table[delta[n], {n, 2, 7}]


I am literally a rookie in the Mathematica software and tbh I have searched online a lot in order to write this piece of code above. Anyway, I am getting an error of the form

$RecursionLimit::reclim: Recursion depth of 1024 exceeded. >>  once I try to calculate the various values of the$μ\$-cycles. (for example: mu[3]).

I would really appreciate your assistance. Thanks!

• Welcome! I suggest the following: 1) As you receive help, try to give it too, by answering questions in your area of expertise. 2) Take the tour and check the faqs! 3) When you see good questions and answers, vote them up by clicking the gray triangles, because the credibility of the system is based on the reputation gained by users sharing their knowledge. Remember to accept the answer, if any, that solves your problem, by clicking the checkmark sign! – user9660 Apr 3 '16 at 17:08
• Seems like you are putting "x=0.5" as your starting point, right? – MathX Apr 3 '16 at 18:08
• @MathX That is correct yes. – Spy93 Apr 3 '16 at 18:13
• your function for mu[n] is recursive, but you didnt define the base case mu[1] anywhere. – AccidentalFourierTransform Apr 3 '16 at 18:36
• @Spy93 I cannot help you further because I dont know anything about the Feigenbaum Constants. About your code I can say: the definition of mu[n] includes mu[n-1]. This means that if you ask Mathematica to calculate mu[10], the result will depend on mu[9] which Mathematica doesnt know its value. Therefore, it will compute mu[9] to get the value, but this depends on mu[8] which Mathematica doesnt know its value - ad infinitum. If at the begining of your code you include mu[1]=1 (for example), then the process ends when the recursion hits the base case mu[1] (hope this helps) – AccidentalFourierTransform Apr 3 '16 at 18:57

Several changes are required to obtain the desired results. First, the syntax error mu[n] == mu must be replaced by mu[n] = mu. Next, initial conditions must be provided for the recurrence in n:

mu[0] := mu /. FindRoot[f[1/2, mu, 2^0] == 1/2, {mu, 0.9}]
mu[1] := mu /. FindRoot[f[1/2, mu, 2^1] == 1/2, {mu, 0.9}]
mu[2] := mu /. FindRoot[f[1/2, mu, 2^2] == 1/2, {mu, 0.9}]


However, even with these changes, FindRoot in lines 4 - 6 of the code in the question does not evaluate mu[n - 1], etc in the initial guesses for its search. This issue can be circumvented by moving the evaluation of these quantities outside FindRoot. Altogether, the revised code becomes

f[x_Real, mu_Real, n_Real] := Block[{y = x}, Do[y = 4 mu y (1 - y), n]; y]
mu[0] := mu /. FindRoot[f[1/2, mu, 2^0] == 1/2, {mu, 0.9}]
mu[1] := mu /. FindRoot[f[1/2, mu, 2^1] == 1/2, {mu, 0.9}]
mu[2] := mu /. FindRoot[f[1/2, mu, 2^2] == 1/2, {mu, 0.9}]
delta[n_] := (mu[n - 1] - mu[n - 2])/(mu[n] - mu[n - 1])
mu[n_] := mu[n] = Module[{mu1 = mu[n - 1] + (mu[n - 1] - mu[n - 2])/delta[n - 1],
mu2 = mu[n - 1] + (mu[n - 1] - mu[n - 2])/(delta[n - 1] + .01)},
mu /. FindRoot[f[1/2, mu, 2^n] == 1/2,  {mu, mu1, mu2}, AccuracyGoal -> 10]]
Table[mu[n], {n, 0, 7}]
Table[delta[n], {n, 2, 7}]

(* {0.5, 0.809017, 0.87464, 0.88866, 0.891667, 0.892311, 0.892449, 0.892478} *)
(* {4.70894, 4.68077, 4.66296, 4.6684, 4.66895, 4.66916} *)

• How does f[1/2, mu 2^n] even evaluate, since 1/2 and 2^n aren't Real? When I run f[1/2, 0.7, 2] it doesn't evaluate, yet these FindRoots work. – Chris K Apr 4 '16 at 2:27
• @ChrisK Interesting observation. Perhaps, FindRoot passes real numbers. In any case, the Real conditions can be eliminated, if Evaluated -> False is added to the FindRoot call. – bbgodfrey Apr 4 '16 at 4:10
• @bbgodfrey Thank you very much for your answer. It was really helpful! – Spy93 Apr 5 '16 at 19:53