# Plotting a table with many y values

I'm trying to plot a bifurcation map for a Lorenz System. I have a function (called estables in the code that follows) which, for a certain value r, returns the list of stationary points for the system with that parameter r.

Then I do this:

r0 = 0;
rf = 10;
dr = 0.5;
ys = Table[estables[ra], {ra, r0, rf, dr}]

ListPlot[ Flatten[ys, {2}], PlotRange -> All , DataRange -> {r0, rf},
PlotStyle -> {Blue, PointSize[Small]}, Joined -> False ]


And I get a weird thing where the points are plotted out of order. That is, the points in the first line are not the ones plotted above the 0, as it should be.

What would be the best way to achieve this?

Edit: This are the differential equations:

ecf = σ (p[t] - f[t]);
ecp = -p[t] + d[t] f[t];
ecd = b (r - d[t] - f[t] p[t]);


and this is the function

t0 = 100;
tf = 105;
dt = 0.05;

estables[ra_] :=

Module[{par, solnum, fs},
par = {σ -> 3., b -> 1, r -> ra};
solnum =
NDSolve[{Derivative[1][f][t] == ecf,
Derivative[1][p][t] == ecp,
Derivative[1][d][t] == ecd,
f[0] == 0.00001, p[0] == 0.,
d[0] == 1} /. par,
{f, p, d}, {t, t0, tf},
MaxSteps -> 100000
];
fs = Flatten[Table[Evaluate[f[t] /. solnum], {t, t0, tf, dt}]];
Union[FindPeaks[fs], -FindPeaks[-fs]]][[;; , 2]]

• Including the definition of estables would certainly help us to help you. May 21, 2016 at 13:08
• I didn't include it for brevity, but yes, you're probably right. May 21, 2016 at 13:12
• Where are t0, tf, and dt? May 21, 2016 at 13:19
• I have a suggestion: why not use WhenEvent[] instead to capture the points where f'[t] == 0? May 21, 2016 at 13:50
• Mmm, I suspected there would be a more sophisticated way to capture the stable points, that looks good. I'm not sure how to use it, though. Which action should be triggered when f'[t] ==0? May 21, 2016 at 13:56

Thanks to J.M. I found a more elegant way to obtain the statinoary values that returned far simpler data, which was easy to plot. This is how my code looks now:

ecf = s (p[t] - f[t]);
ecp = -p[t] + d[t] f[t];
ecd = b (r - d[t] - f[t] p[t]);

par = {s -> 3., b -> 1, r -> 100};
t0 = 800;
tf = 850;

estables[ra_] :=

Module[{par, solnum,  points},
par = {s -> 3., b -> 1, r -> ra};
{solnum, points} = Reap@NDSolve[{
Derivative[1][f][t] == ecf,
Derivative[1][p][t] == ecp,
Derivative[1][d][t] == ecd,
f[0] == 0.00001, p[0] == 0.,
d[0] == 1,
WhenEvent[{f'[t] == 0 && t > t0}, Sow[f[t]]]}
/. par, {f, p, d}, {t, t0, tf}, MaxSteps -> 100000];
Flatten[points]
];

r0 = 0;
rf = 105;
dr = 0.2;

ys = Table[{ra, e}, {ra, r0, rf, dr}, {e, estables[ra]}];
ListPlot[ Flatten[ys, 1], PlotRange -> All , DataRange -> {r0, rf},
PlotStyle -> {Blue, PointSize[0.002]}, Joined -> False,
AxesLabel -> {r, f} ]

• Great job on self answer! When I copy and run your code it appears that points 2 and 3 generate an error. The rest of the points appear to be OK. May 22, 2016 at 13:50
• Actually, I've just found out that my code is wrong. It gives a close enough diagram, but a lot of "spurious" points appear that shouldn't be there. Is that what you mean? I've almost managed to fix it, I will edit in a moment. May 22, 2016 at 20:38