Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

This is what I am stuck on

f[a_][x_] := x^2 + a;
w[a_] := Take[NestList[f[-1.8], 0, 1000], 100]

Now I want to evaluate this for different values of the constant $a$ from -1.8, to 0.3, store the results in a table, suppress the results because they will be very long, and then make a plot of $a$ versus $x$ using ListPlot.

How can I evaluate for different values of $a$ and then plot?

I have tried the Do function but can't get it to work.

share|improve this question
    
Are you really expecting overflows to occur here? –  J. M. May 4 '12 at 6:06
2  
@JohnStamos, could you perhaps make your title a little more informative? –  user21 May 4 '12 at 6:08
1  
There is no $c$ in your code. Essentially what you are trying to do is work out the Mandlebrot set just along the real line, which is the boring part of the set. It might help to have some more explanation of what you are trying to do. –  Verbeia May 4 '12 at 6:36

2 Answers 2

up vote 3 down vote accepted

There is no a in your definition of the function w[], so it's no wonder that you can't get it to work.

You could try the following, but you would hit the same Overflow problems as in your previous question.

w[a_] := Take[NestList[f[a], 0, 1000], 100]

I would suggest an alternative approach using NestWhile:

w[a_] := NestWhile[f[a], 0, Abs[#] < 1000. &, 1, 100]

You will only get the "last" element in the iteration, defined either as the 100th iteration - it will have converged by then if it was going to - or the last iteration before the absolute value of the result broke through 1000 (i.e. it was diverging). That way, you break out before the iteration overflows.

You can then get your list of results using the Table function (notice I've also preserved the values of $a$, which I've labeled $i$ for the purpose of defining the iteration in the table, so the values of the $a$ parameter are shown correctly in the subsequent plot:

result=Table[{i, w[i]}, {i, -0.2, 1.8, 0.01}];

And plot using ListLinePlot:

ListLinePlot[result, PlotRange -> All, Frame -> True, AxesOrigin -> {0, 0}]

enter image description here

One thing I would mention is that it's not clear whether you should focus on a single starting point. You could instead define:

ff[a_, x0_] := NestWhile[#^2 + a &, x0, Abs[#] < 100. &, 1, 100]

And construct your results as:

result = Table[{i, ff[i, 0.1]}, {i, -0.2, 1.8, 0.01}];

You could then even have a two dimensional table:

result2 = Flatten[
  Table[{i, j, ff[i, j]}, {i, -0.2, 1.8, 0.01}, {j, -0.2, 0.9, 0.1}],  1];

Which you could visualise using ListDensityPlot:

ListDensityPlot[result2]

enter image description here

share|improve this answer
    
Thats a pretty awesome reply! thanks –  John Stamos May 4 '12 at 7:06
    
@JohnStamos you are welcome. As described in the FAQ, it is probably best to wait to see if a better answer is posted before accepting mine. Since you're new it would be a good idea to familiarise yourself with how StackExchange works, including the editing, voting and accepting features. And welcome to the Mathematica site. –  Verbeia May 4 '12 at 7:18
    
when i try to plot it though i get the error "result is not a list of numbers"or a pair of numbers" –  John Stamos May 4 '12 at 7:21

It sounds to me like you're trying to plot a bifurcation diagram.

f[a_][x_] := x^2 + a;
eventualOrbit[a_] := 
  Drop[NestList[f[a], 0.0, 500], 100];
column[a_] := {a, #} & /@ eventualOrbit[a];
bifPic = Graphics[{
   PointSize[0.001], Opacity[0.1],
     Point[Flatten[Table[column[a], {a, -2, 0, 0.001}], 1]]
   }, AspectRatio -> 1/2]

enter image description here

We needn't worry about divergence here, since the orbit of zero happens to be bounded for all $a$ in this region. There are curves evident in the image, called the critical curves defined by $Q_n(a) = f_a^n(0)$. We can plot them on top of the bifurcation diagram like so:

F[x_]  = Nest[f[a], x, 4];
a4 = a /. FindRoot[{F[x] == x, F'[x] == -1}, {x, .8}, {a, -1.4}];
Q[n_, a_] := Nest[f[a], 0, n];
curves = Plot[Table[Q[n, a], {n, 1, 7}],
   {a, -2, a4}, PlotPoints -> 100, 
   PlotStyle -> {Thickness[0.0025], GrayLevel[0.2]}];
Show[{bifPic, curves}]

enter image description here

share|improve this answer

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.