I encountered a definite integral (with parameter $k$) that produces the "chaotic" figures below.
$$F(k) \equiv \int_{\theta_k}^{\theta_{k+1}} f(\theta, J \sin\theta - k)\,\mathrm{d}\theta~, \qquad f(\theta,u) \equiv \left( \frac{2u - L\sin\theta }{1 - L \sin\theta} \right)^2 + \left( \frac{2 - 2u + L \sin\theta }{1 - L \sin\theta} \right)^2\\ \theta_k \equiv \sin^{-1}\bigl( \frac{k}{J} \bigr)~, \qquad |L|<1,\quad L\neq0,\quad J>k+1 \geq1,\quad J \in \mathbb{R},~k \in \mathbb{N}$$
It is a line integral along a sinusoidal path $u = J \sin\theta - k$, going from $u = 0$ to $u = 1$, where the angles are $\theta_k$ to $\theta_{k+1}$.
My Question:
Is this a bonafide bifurcation to chaos? Is it some other mathematical phenomenon? Or is it likely an implementation flaw?
Each of the plots is $F(k)$ versus $k$ for a distinct large $J$. Basically my intention was to use $\frac{k}J$ to approximate 0 to 1.
Roughly speaking, the discrete points start out nicely following a curve and then spreads out for perhaps $J > 10^4$.
There is a parameter $L$ which the distribution of the points depends on.
The actual parameters of these 4 plots are not important cuz this bifurcating-like behavior is very common. Following the figures is the code block of a Manipulate
module that one can play with and generate similar plots
With[{Mn = 0.01, Mx = 16, SpL = Spacer@15, SpS = Spacer@3, Fz = FieldSize -> 2.5},
DynamicModule[{Opt, indef, Jp = 5(* Log@J *), Ns = 3(* Log@ sample points *), ShowAll = False, ImgSz = 200, PtSz = 1/100},
indef[t_, J_, L_, k_] := 1/(4 \[Pi]) (3 (1 - L^2)^(3/2))/(1 + 3 (1 - L^2)^(3/2))Divide[ 4/Sqrt[1 - L^2] ( (2 J - L (1 + 2 k)) ((2 J - L) (3 L^2 - 2) - 2 L^3 k) ArcTan[(L - Tan[t/2])/Sqrt[1 - L^2]] + (2 J + L (1 + 2 k)) (2 L + 2 J (3 L^2 - 2) + L^3 (2 k - 1)) ArcTan[(L + Tan[t/2])/Sqrt[1 - L^2]]) + 1/(-1 + L^2 Sin[t]^2) (-4 (2 J - L) L (2 - 3 L^2 + L^4) (1 + 2 k) t - L (-4 J L (4 - 5 L^2 + L^4) + 4 J^2 (8 - 5 L^2 + L^4) + L^2 (-5 L^2 + L^4 + 8 (1 + 2 k (1 + k)))) Cos[t] + L^3 ((-2 J + L) (-1 + L^2) (-4 (1 + 2 k) t Cos[2 t] + (-2 J + L) Cos[3 t]) + 8 J (1 + 2 k) Sin[2 t]))
, 2 L^3 (-1 + L^2)];
Manipulate[ (* lazy foolproofing *) If[L == 0, L = Mn (-1)^RandomInteger[]];
Opt = {PlotStyle -> PointSize -> Dynamic@PtSz, ImageSize -> Dynamic@ImgSz, Joined -> False, Filling -> False};
Row@{ DiscretePlot[{ indef[ArcSin@((k + 1)/J), J, L, k] - indef[ArcSin@(k/J), J, L, k], ArcSin@((k + 1)/J) - ArcSin@(k/J) }, {k, 0, -1 + Floor@J, Max[1, Floor[J]/10^Ns]}, Evaluate@Opt, PlotRange -> If[ShowAll, All, Automatic]],
DiscretePlot[ Divide[indef[ArcSin@((k + 1)/J), J, L, k] - indef[ArcSin@(k/J), J, L, k]
, ArcSin@((k + 1)/J) - ArcSin@(k/J) ], {k, 0, -1 + Floor@J, Max[1, Floor[J]/10^Ns]}, Evaluate@Opt]}
, Column@{ Row@{"Max Log@J", SpS, InputField[Dynamic@Jp, Fz], SpL,
"ShowAll", SpS, Checkbox@Dynamic@ShowAll}, Row@{"Log@Ns", SpS, SetterBar[Dynamic@Ns, Range@4]}, Row@{"ImgSz", SpS, InputField[Dynamic@ImgSz, Fz], SpL, "PtSz", SpS, InputField[Dynamic@PtSz, Fz]} }
, {{J, 29000}, Dynamic[10^(Jp - 2)], Dynamic[10^Jp],
Dynamic[10^(Jp - 2)]}, {{L, -0.17}, -1 + Mn, 1 - Mn, Mn}, TrackedSymbols :> {Jp, J, L, Ns, ShowAll, PtSz}, ControlPlacement -> Left] ] ]
The codes above also generates a plot comparing the blue curve $F(k)$ with $\theta_{k+1}-\theta_k$ in yellow.
How $F(k)$ is implemented in the module above: hard-code the indefinite integral $\int f$, which is previously computed outside, then evaluate it at $\theta_{k+1}$ and $\theta_k$ then take the difference.
The code block below is how I compute the indefinite integral and FullSimplify
to improve it slightly.
ClearAll[k, J, u, f, indef, indefSimp, L, myAss];
f[t_, u_, L_] :=((2 u - Jc Sin@t)/(1 - Jc Sin@t))^2 + ((2 - 2 u + Jc Sin@t)/(1 + Jc Sin@t))^2;
myAss = -1 < L < 1 < k <= J && K != 0; Timing[indef = Integrate[Sin@t f[t, J Sin@t - k, Jc], t, Assumptions -> myAss]]
Timing[indefSimp = FullSimplify[indef, myAss && 0 < t < Pi/2]] (* this takes about 30 seconds up to 2 minutes *)
Just in case if this is relevant, I also have a Plot3D
module showing the surface of $f(\theta,u)$ like these figures below (left: $L = -0.86$ and right: $L = 0.23$ arbitrary).
Note that $f(\theta,u)$ is a legitimate probability density over $\theta \in [0,2\pi)$ and $u \in [0,1)$ with a normalization constant that is actually implemented in the DiscretePlot
Manipulate
above (though irrelevant I think there's a multiplicative factor of 4 missing).
With[{ mM = Max[ #[[1]], Min[ #[[2]], Re@#[[3]] ] ] &, SpL = Spacer@15, SpS = Spacer@3 (* Spacer Large & Small *), Mn = 0.01, Mx = 10, Fz = FieldSize -> 2.5},
DynamicModule[{ImgSz = 250, Opt, intP = True}, Manipulate[{L, p} =
mM /@ {{-1 + Mn, 1 - Mn, L}, {If[intP, 1, 0], Mx, If[intP, p, Ceiling@p]}}; Opt = {BoxRatios -> {Pi, 2, 2}, PlotPoints -> {4, 1} 10,
ClippingStyle -> None, ImageSize -> Dynamic@ImgSz};
Plot3D[ ((2 u - L Sin@t)/(1 - L Sin@t))^p + ((2 - 2 u + L Sin@t)/(1 + L Sin@t))^p , {t, 0, 2 Pi}, {u, 0, 1}, Evaluate@Opt],
Row@{"ImgSz", SpS, InputField[Dynamic@ImgSz, Fz], SpL, "Integer p", SpS, Checkbox@Dynamic@intP}, {{L, 0.23}, -1 + Mn, 1 - Mn, Mn}, {{p, 2}, Dynamic@If[intP, 1, 0], Mx, Dynamic@If[intP, 1, Mn]}, TrackedSymbols :> {L, p, intP}, ControlPlacement -> Left] ] ]
This module of Plot3D
is an overkill that covers the $f(\theta,u)$ in question, as a special case where the power is $p=2$, as well as general $p$.
Optional follow up Non-Mathematica${}^{\text{TM}}$ mathematical question
If this is a real mathematical behavior, then it seems that for example $\lim_{J \to \infty} F(k)/(\theta_{k+1}-\theta_k)$ cannot be taken. What would be the proper analysis to bound $F(k)$ for large $J$?
I actually have seen something like this before, but honestly I'm not sure if mine is an example of the ubiquity of chaos or if it's an artifact. Originally I was dealing with something that should converge (the blue dots should be a smooth curve). I cannot make the call if there's something overlooked in my analytic work or is there some unintentional effects due to how I code. For example, $F(k)$ should be non-negative but in the plots there's about half below zero. Any input will be appreciated.
Please let me know if this post can be improved (too few info, too much info, etc). Thank you for your time.
WorkingPrecision
, so you could change that value and see if it affects the plots? $\endgroup$WorkingPrecision
in my plot option ofDiscretePlot
, varying from 6 to 1000. It doesn't appear to affect the outcome much, if at all. Maybe I'm doing this wrong. I'll test some more. $\endgroup$