# Bifurcation Diagram for x' = (rx) - sinx

I solved for finding the fixed points analytically for various ranges of u for x'[t] = (ux) - sinx

1. for u = 0
2. 0 <u <<1
3. u>1
4. 0>u>-infinity

and its a very complicated bifurcation diagram as shown in this picture below. The dotted lines represent unstable region and stable lines are shown by the normal lines.

I am trying to plot it, however only empty plots get returned, and it gives me error for using NDSolve. I am also a beginner in using mathematica. Any help would be appreciated.

    Clear["Global*"]

dxdt[x_, u_] := u*x - Sin[x]

xInitial = 1;

xSolution = NDSolveValue[{x'[t] == dxdt[x[t], u], x[0] == xInitial}, x, {t, 0, 1},
Method -> {"StiffnessSwitching", "NonstiffTest" -> Automatic}];

stablePoints = {};
unstablePoints = {};

Do[
sol = NDSolveValue[{x'[t] == (D[dxdt[x[t], u], u] /. x -> xSolution), x[0] == xInitial},
x, {t, 0, 1}];
If[Chop[sol[1]] == 0,  (* Chop is used to handle numerical inaccuracies *)
AppendTo[stablePoints, {u, xSolution[u]}],
AppendTo[unstablePoints, {u, xSolution[u]}]
],
{u, -2, 2, 0.01}
];

Show[
ListPlot[stablePoints, PlotStyle -> Black],
ListPlot[unstablePoints, PlotStyle -> Black, PlotMarkers -> {"*", 10}],
Plot[xSolution[u], {u, -2, 2}, PlotRange -> All,
AxesLabel -> {"u", "x"},
PlotLabel -> "Bifurcation Diagram for x' = u*x - Sin[x]"]
]


And its outputs

• NDSolve can not deal with symbolic parameters. Therefore you get an error because of x'[t] == dxdt[x[t], u] where u is symbolic. Commented Mar 12 at 9:23
• Are you trying to get the plot, ContourPlot[0 == u x - Sin[x], {x, -6 Pi, 6 Pi}, {u, -1.1, 1.1}, MaxRecursion -> 3]? Alternatively, ParametricPlot[{Sin[x]/x, x}, {x, -6 Pi, 6 Pi}, AspectRatio -> 1, PlotRange -> All]? Commented Mar 12 at 13:48
• @Goofy, I am getting a plot using ContourPlot but it does not differentiate between unstable and stable fixed points by representing unstable using dotted lines. I tried to incorporate that but it outputs everything belonging to the stable subspace. I also printed the stable points and unstable points separately as (x, u) but their values of the tuples are extremely close, so it becomes continuous. Commented Mar 12 at 15:04
• Do you mean like this?: ParametricPlot[{Sin[x]/x, x}, {x, -6 Pi, 6 Pi}, MeshFunctions -> {Function[{u, x}, x Cos[x] - Sin[x]]}, MeshShading -> {Dotted, Automatic}, Mesh -> {{0}}, AspectRatio -> 1, PlotRange -> All] Commented Mar 12 at 21:21
• @Moo It's easy to fix, no? Why don't you do the work for them? Or you could fix their approach. Before perfecting my approach, one might consider they are asking about a different approach, which is probably more robust (and complicated), and that my ContourPlot[] approach is probably better than ParametricPlot[], although MeshShading doesn't work for it. Commented Mar 14 at 15:45

This might not be the plot you're after, but I think it shows what you want to show:

myArrowheads[c1_, c2_] :=
Table[{s, Automatic, Graphics@N@{c1, EdgeForm[c2], Disk[]}},
{s, ConstantArray[0.005, 5]}]];
splitLine = GraphicsComplex[c_, g_, o___] :> GraphicsComplex[c,
g /. Line[p_] /; Max@c[[p, 1]] > 1.1 :>
{Line /@ GatherBy[p, c[[#, 1]] - 1 > 0 &]}, o];
Show[
ContourPlot[
ConditionalExpression[u  x - Sin[x], u - Cos[x] <= 0] == 0
, {u, -0.7, 1.5}, {x, -6 Pi, 6 Pi}
, MaxRecursion -> 4] /.
ColorData[97][1]], Arrow@pts},
ContourPlot[
ConditionalExpression[u  x - Sin[x], u - Cos[x] >= 0] == 0
, {u, -0.7, 1.5}, {x, -6 Pi, 6 Pi}
, MaxRecursion -> 4, ContourStyle -> Dotted] /. splitLine /.

As usual, good plots are a lot of work. Better spacing between the dots, for instance. One might extract the lines and apply the approaches in How can I resample a list of {x,y} data for evenly spaced points? or Divide a Graphics Line into $i$ equal lengths One needs to pay attention to the distance in the actual geometry versus the distance in the scaled screen/print image with different horizontal and vertical scales. Also, if you need the points where u == 0, you will have to figure how how to space the other dots around them. Or if you want only those points, then you can add them with Epilog or Show[plot, Graphics[<code for dots>]]`.