8
$\begingroup$

So I have the general code from the PoincareSection documentation that is changed up for a Driven Damped Pendulum:

coupledDiffEq = 
 {ω'[t] == -(1/q) ω[t] - Sin[θ[t]] + 
  g*Cos[ϕ[t]],
  θ'[t] == ω[t],
  ϕ'[t] ==  Drive};

data = Block[{q = 3.9, g = 1.5, Drive = 1}, 
Reap[NDSolve[{coupledDiffEq, 
θ[0] == 0, ω[0] == 0, ϕ[1] == 2*Pi, 
WhenEvent[Mod[ϕ[t], (2*Pi)] == 0, 
Sow[{θ[t], ω[t]}]]}, {}, {t, 0, 100000}, 
MaxSteps -> ∞]]][[-1, 1]];
ListLinePlot[data, PlotRange -> {{-3, 3}, All}]

Then I ListPlot data and it comes out as not a PoincareSection. Often it looks like a phase plot more than anything else. I have attempted to change the

Mod[ϕ[t], DriveFrequency]

where DriveFrequency should be Drive / (2*Pi), but that often does nothing. The system of equations evaluates and produces points, they just are random though. How do I get a nice Poincaré section with this?

my result

I am getting this when I should be getting

correct Poincaré section

$\endgroup$
1
  • $\begingroup$ Yes, when I use the Parameters (q = 3.9, g =1.5, Drive = .6667) I get something that looks more like a phase diagram. Imager Link to result: imgur.com/8natnP0 $\endgroup$ – Escap3faith Nov 7 '15 at 23:12
8
$\begingroup$

Update

J.M. has some good comments, so I'm going to implement them.

coupledDiffEq = {θ''[t] == -(1/q) θ'[t] - 
 Sin[θ[t]] + g*Cos[ϕ[t]], ϕ'[t] == drive};

data = Block[{q = 3.9, g = 1.5, drive = 2/3}
   , First@Last@Reap@NDSolve[{
        coupledDiffEq, θ[0] == 0, θ'[0] == 0, ϕ[1] == 2 π
        , WhenEvent[Mod[ϕ[t], 2 π] == 0, Sow[{Mod[θ[t], 2 π, -π], θ'[t]}]]
        }, {}, {t, 0, 100000}, MaxSteps -> ∞]
   ];
ListPlot[Drop[data, 100]]

I dropped the first 100 points to make sure that we aren't capturing any of the transients.


I played around with your code, and couldn't get it to do what I wanted, so I went this route instead. Note that you need to move θ[t] back into the interval to π in order to get the attractor, and you should not use ListLinePlot, because you will get a mess.

ClearAll[coupledDiffEq, ω, q, θ, g, ϕ, drive]
coupledDiffEq := {ω'[t] == -(1/q) ω[t] - Sin[θ[t]] + g*Cos[ϕ[t]], θ'[t] == ω[t], ϕ'[t] == drive};

sols = Block[{q = 3.9, g = 1.5, drive = 2/3}
  , First@NDSolve[
     {coupledDiffEq, θ[0] == 0, ω[0] == 0, ϕ[1] == 2*π}
     , {θ[t], ω[t]}
     , {t, 0, 100000}
     , MaxSteps -> ∞
    ]
 ];

 data = Table[{Mod[θ[t] + π, 2 π] - π, ω[t]} /. sols, {t, 20*(2 π)/(2/3), 100000, (2 π)/(2/3)}];
 ListPlot[data, PlotRange -> All]

enter image description here


Update

This code should work too, more along the lines of what you did:

coupledDiffEq = {ω'[t] == -(1/q) ω[t] - Sin[θ[t]] + g*Cos[ϕ[t]], θ'[t] == ω[t], ϕ'[t] == drive};

sols = Block[{q = 3.9, g = 1.5, drive = 2/3}
  , First@Last@Reap@NDSolve[{coupledDiffEq, θ[0] == 0, ω[0] == 0, ϕ[1] == 2 π
  , WhenEvent[Mod[ϕ[t], 2 π] == 0, Sow[{Mod[θ[t] + π, 2 π] - π, ω[t]}]]
 }, {θ[t], ω[t]}, {t, 0, 100000}, MaxSteps -> ∞]];

ListPlot[sols]

enter image description here

$\endgroup$
6
  • 4
    $\begingroup$ Two notes: 1. As I noted in this answer, you can pass an empty list as the second argument of NDSolve[] so that it doesn't save the InterpolatingFunction[] objects, potentially saving memory. OP seems to have tried this, but could not make it work. 2. Mod[θ[t] + π, 2 π] - π is slightly more compactly done as Mod[θ[t], 2 π, -π]. $\endgroup$ – J. M.'s torpor Nov 8 '15 at 8:32
  • $\begingroup$ You are amazing! Its quite interesting how you have to play around with Sow[...] and θ[t] to get this to work! I really appreciate your help! $\endgroup$ – Escap3faith Nov 8 '15 at 14:20
  • 1
    $\begingroup$ @Escap3faith, BTW, you don't really need ω for your system (assuming you only used it as a temporary intermediate); NDSolve[] can handle θ''[t] == -θ'[t]/q - Sin[θ[t]] + g Cos[ϕ[t]] perfectly, and your Sow[] becomes Sow[{Mod[θ[t], 2 π, -π], θ'[t]}]. $\endgroup$ – J. M.'s torpor Nov 8 '15 at 17:28
  • $\begingroup$ @J.M. Thanks for comments and links! I actually couldn't think of a good reason why not to actually solve for the functions, but not returning a memory-intensive InterpolatingFunction is a good one. I will modify the post with your comments, I think, to make the answer more worthwhile. $\endgroup$ – march Nov 8 '15 at 17:38
  • 1
    $\begingroup$ @J.M. I don't think I do, but I was reading a paper and they broke it down into 3 equations due to the necessity of ϕ[t], and I just followed out of frustration for the project. I'm 100% positive I don't need all 3 coupled equations because Mathematica has the power to deal with them, but rather just the 2. Paper: thphys.uni-heidelberg.de/~gasenzer/… $\endgroup$ – Escap3faith Nov 8 '15 at 19:50

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.