Table of Contents
perorb - interactive/iterative search for XY-periodic orbits in symmetric
potentials
perorb out=orbit [parameter=value]
perorb
searches for periodic orbits in two ways:
In the iterative mode it will
integrate an orbit for a half or whole period, and iteratively converge
towards a periodic orbit assuming regular phase space. It uses a convergence
criterion in the surface of section coordinates. It then can proceed, and
attempt to find family members with slightly different energies or spatial
extent. The orbits can be stored in a binary file along the way, and also
a table with major an minor axis rotation curves, period and energy can
be stored in a table.
In the interactive mode the selected orbits are integrated
for a specified number of surface of section (S.O.S.) crossings, and their
S.O.S. coordinates are plotted; from a visual inspection the central periodic
orbit(s) can then be estimated at this energy, and used to get good estimates
in the iterative mode to find the periodic orbits.
This technique has been
used in various forms by various researchers, and is described in its
most original form in Henon (1965, Ann.Astr. 28, 499). See also El-Sabaa &
Sherief (1990) Ap&SS 167, 305 for some more recent work.
The following
parameters are recognized in any order if the keyword is also given:
- in=in-file
- An optional input file in standard orbit(5NEMO)
format. If given, the input
conditions are taken from the first phase point in this orbit (using the
keyword last= you can select the last orbit). Default: none. out=out-file
output file, into which the orbits are accumulated. The file is in standard
orbit(5NEMO)
format. If no filename is given, the iterative mode is followed.
[Default: none].
- freqout=rate
- Frequency for output of the timeslots for
the periodic orbits. Only used in iterative mode. [default: 1].
- maxsteps=number
- maximum number of reserved integration steps. [default: 5000].
- dt=time-step
- time step (constant) to be taken [default: 0.05].
- phase=trial-phase
- Initial
trial phase space coodinates for orbit. Although currently only planar orbits
in the XY-plane can be investigated, all 6 phase-space coordinates need
to be specified. There should be only two non-zero entries in this case as
the launch coordinates. If a 7th number is given (optional), this is taken
as the trial energy for the first orbit. The trial launch velocity is then
calculated from the available potential at the specified launch position.
- step=step,perturb
- Step in phase space to take to get next orbit (the step
is taken in the direction given by the dir keyword. The second number is
the perturbation in direction dir used to start the first iteration that
should result in finding the periodic orbit. Default: 0.1,0.01
- dir=x|y
- Direction
in which to step for next orbit. This also determines in which surface of
section the searching for the periodic orbit is done. Only x and y are allowed
for now, since only planar XY-orbits are handled properly.
- norbit=n
- Number
of orbits to investigate. [Default: 1]
- accuracy=eps
- A relative accuracy
parameter to which the final surface of section coordinate is to be accurate.
Only used in iterative mode. [Default: 0.001].
- maxiter=
- Maximum iterations
allowed to find a periodic orbit. Default: 50.
- ncross=count
- Number of crossings
through surface of section used to plot/iterate. This count does include
the half-period crossings if the keyword period=2 is choosen. [default:
1]
- period=2|1
- The number of surface of section coordinates to generate per
period. Must be either 1 or 2. A value of 2 make search faster, but a bit
more difficult for irregular orbits or near resonances. Also non-symmetric
orbits cannot be found this way. [default: 2].
- potname=name
- name of file
of potential(5NEMO)
descriptor [default: plummer].
- potpars=par-list
- List
of parameters to the potential descriptor. The first parameter must be the
pattern speed in the x-y plane, although rotating frames of reference are
not yet supported. The remaining parameters are used by the _inipotential()
routine in the potential descriptor. [default: none - let them be defined
by routine itself].
- potfile=name
- name of an optional datafile to the potential
descriptor. This might be an N-body snapshot or list of spline fit coefficients
etc. [default: none].
- tab=table_file
- If set to a filename, a table, which
normally appears on stdout is now written to a file. It contains the two
axes rotation curves, number of steps in the orbit, numer of iterations
in surface of section, the full period and the energy: pos,vel,pos,vel,nsteps,niter,period,energy.
[Default: none].
- mode=int_mode
- Specify the integration mode. Any one of euler,
leapfrog, rk2 or rk4 can be given. [Default: rk4].
- last=t|f
- Controls which
orbit from the input file is used for the initial conditions. By default,
the first one, otherwise the last. Default: f.
- headline=
- Random user verbiage
added to the output file. Default:none
The potential (potname= etc)
needs to have a proper potential defined for this routine to work.
By
default, perorb will find a circular orbit in a plummer potential near
R=1:
perorb
INIPOTENTIAL Plummer: [3d version]
### Warning [perorb]: DRY RUN: No output orbit file created
# x0 v0 y1 x1 NPT NITER PERIOD ETOT LZ_MEAN
ETOT_ERR
1.007320 0.595679 1.007323 -0.595674 533 3 10.625226 -0.527107 -0.600039 4.08912e-06
The printed columns are pos1(t=0), vel2(t=0), pos2(t=P/2), vel1(t=P/2),
l, iter, Period, Energy, meanAngMom. See also orbstat(1NEMO)
for another
view of the statistics of these orbits, showing the curious property of
the 1:1 resonance orbit in a logarithmic potential
% perorb potname=log potpars=0,1,1,0.9 norbit=10 out=log.orb
# x0 v0 y1 x1 NPT NITER PERIOD ETOT LZ_MEAN
ETOT_ERR
0.956311 0.675829 0.629557 -1.140210 280 5 5.560441 0.984048 -0.679534 -1.4864e-05
1.020431 0.731892 0.706384 -1.174306 290 5 5.761082 1.094431 -0.785260 -1.55914e-06
1.086113 0.781263 0.781720 -1.205717 300 5 5.972464 1.203861 -0.892278 -1.50248e-05
1.153376 0.825038 0.856210 -1.234575 311 5 6.194336 1.312000 -1.000653 -7.58749e-06
1.222193 0.864064 0.930308 -1.260978 323 4 6.426479 1.418565 -1.110458 -1.37877e-05
1.292542 0.898976 1.004360 -1.285151 335 4 6.668597 1.523326 -1.221855 -7.93163e-06
1.364403 0.930272 1.078516 -1.307186 348 4 6.920329 1.626100 -1.334623 -1.17395e-05
1.437707 0.958411 1.153092 -1.327410 361 4 7.181385 1.726745 -1.448891 -8.47802e-06
1.512414 0.983745 1.228047 -1.345782 374 4 7.451368 1.825161 -1.564586 -9.25693e-06
1.588487 1.006573 1.303499 -1.362555 388 4 7.729841 1.921275 -1.681411 -8.9171e-06
% orbstat log.orb
# T E x_max y_max u_max v_max j_mean j_sigma
0 0.984048 0.629601 0.956311 0.675829 1.14031 -0.679593 0.0252821
0 1.09443 0.706429 1.02043 0.731892 1.1744 -0.785327 0.0292126
0 1.20386 0.781738 1.08611 0.781263 1.20575 -0.892351 0.0331657
0 1.312 0.85621 1.15338 0.825038 1.23457 -1.00073 0.037188
0 1.41857 0.930308 1.22219 0.864064 1.26098 -1.11054 0.0412889
0 1.52333 1.00436 1.29254 0.898976 1.28515 -1.22194 0.0454228
0 1.6261 1.07857 1.3644 0.930272 1.30726 -1.33472 0.0496375
0 1.72675 1.15309 1.43771 0.958411 1.32741 -1.44899 0.0538812
0 1.82516 1.22807 1.51241 0.983745 1.34581 -1.56469 0.0581483
0 1.92128 1.30352 1.58849 1.00657 1.36259 -1.68152 0.0624934
% orbstat log.orb | tabmath - - '(%4*%5)/(%3*%6)' all
0.900218
0.900214
0.900228
0.900224
0.900222
0.90022
0.900205
0.900233
0.900213
0.900214
The last command shows the ratio of the angular momentum along the major
and minor axis, and that they equal the axis ratio (q) of the potential.
orbstat(1NEMO)
, potlist(1NEMO)
, orbint(1NEMO)
, newton0(1NEMO)
,
potential(5NEMO)
Peter Teuben
pjt/orbit original sources
22-may-90 V1.0 created from old Cyber program 'PERORB' PJT
24-may-91 V1.1 rotating XY-frames, fixed energy option PJT
19-apr-95 V1.5 various, rk4 is now default integrator PJT
1-mar-03 minor code cleanup, less lies in the man page, added etot_err PJT
19-aug-04 V1.6 added last= and fixed an allocation problem PJT
Table of Contents