Re: Poincare section for double pendulum

*To*: mathgroup at smc.vnet.net*Subject*: [mg120596] Re: Poincare section for double pendulum*From*: JUN <noeckel at gmail.com>*Date*: Sat, 30 Jul 2011 06:02:05 -0400 (EDT)*Delivered-to*: l-mathgroup@mail-archive0.wolfram.com*References*: <j0u7jr$fbg$1@smc.vnet.net>

On Jul 29, 5:04 am, gal bevc <gal.b... at gmail.com> wrote: > Hello, > > I'm a relatively new user of Mathematica, who doesn't have much of > programming skills. For my undergraduate assingment I must analyze chaotic > motion of double pendulum. > > Until now i have got system of differential equations for equations of > motion for double pendulum(i have x''[t]=function(t) and > y''[t]=function(t)). System of differential equations can be solved for 4 > inital conditions, x[0],y[0],x'[0] and y'[0]. With using function NDSolvei > got functions of angles and angular velocities for upper and lower pendulum > with respect to time, x[t],x'[t],y[t] and y'[t]. > > To get a poincare section of double pendulum, i have to record position of > y[t] and y'[t] whenever x[t] is equal to zero and the velocity of x'[t] is a > positive number. In the end I must get some sort of phase diagram y[t] and > y'[t]. > Because this is a Hamilton non-dissipative system, inital energy of the > system is a constant of time and initial energy is a function of initial > conditions. To get a real poincare diagram i must repeat the procedure > described above for different initial conditions, but for the same energy > level. I need mathematica to use some random numbers for initial conditions > in a way that the initial energy of the system stays the same. So i must > repeat procedure for poincare section(surface of section) for let's say 50 > different initial conditions and then display all results in one y[t],y'[t] > diagram. > Hope that someone can help me. > > Thank you, > Gal Bevc Hi, one way of doing this is by using StepMonitor during the numerical solution to look for zero-crossings of the variable x[t]: (a) Define an empty list, say zeros = {}; This will be used to collect the zeros in your time interval. (b) Define an auxiliary variable "lastX" and set it equal to the initial value of x (say xInitial), lastX = xInitial; (c) In your NDSolve block, add the following: ...,StepMonitor:>(If[x[t]*lastX<0,AppendTo[zeros,t]];lastX=x[t]) This tests for zero crossings of x[t] between the current time step and the previous one. (d) When NDSolve is done, the list "zeros" contains a set of approximate t values with x[t]=0 that you can then use to refine using FindRoot. Let's say your solution (in the form of rules {x->..., y->...}) is stored in "sol", then say Map[FindRoot[Evaluate[x[t]/.sol],{t,#}]&,zeros] That should give you the values of t at which you would then evaluate y[t] and y'[t] to make the points for the Poincare section. Jens