Rocking block integrator
Download
single_block_integrator.tar.gz
How to execute
Assume you have successfully compiled the source. It has to be
executed from the command line, just typing:
./main_single_block x0 y0 t0 epsilon numiterates maximpacts
where:

(x0,y0,t0) is the initial condition. Note that time is included as variable.

epsilon is the perturbation parameter, the magnitude of the periodic forcing

numiterates (integer) is the maximum number of iterations of the stroboscopic map that we want to compute

maximpacts (integer) is the maximum number of iterations of the impact map that we want to compute
The values of the frequency of the forcing and the restitution
coefficient have to be set in the file "system_params.dat".
What will (should) the program do?
It will integrate the system until one of the following things occur:

The number of iterations by the strobosopic map reaches numiterates

The number of iterations of the impact map reaches maximpacts

The trajectory escapes to infinity. There is an internal parameter
("stop_tol") to decide this (see below).
Note that the first condition, 1., 2. or 3., that occurs will stop the
program. Hence, if one wishes necessary to compute an given number of
iterations of the stroboscopic map then set "maximpacts" to a very
high value to ensure that 1) occurs first (similarly for the impact
map).
What will the program return?
The program will create three files:

stroboscopoicmap.dat: Contains the iterates of the stroboscopic map

trajectory.dat: contains points to plot the timecontinuous trajectory.

impactmap.dat: contains the iterates of the impacts map.
All three files will contain three columns: x y t
Note that these three files are removed at each execution. Make copies of them if you want to keep data.
How can I plot the results?
Gnuplot is an excellent option, but matlab would do the trick to.
There may be thousands of other options.
What parameters can be (easily) changed/tuned?
Basically all of them. These are set in the preamble of the file
"main_singleblock.cpp", and refer to plotting options, integration
options, tolerances for the Newton method and stopping conditions.
Note that if any is changed then to the program has to be recompiled.
All of them are described in the preamble of the main program. Of
course, the frequency (omega) and restitution coefficient (r), which
are given in the system_params.dat file can be changed at any time
without recompiling. Also the perturbation parameter (epsilon), the
initial conditions and the maximum number of iterates for the maps can
be also changed as desired, as they are the input parameters when
calling the main program. Note that, in principle, there is no need
to change anything in the "integrate_singleblock.cpp" file, which
contains the machinery.
Where are the differential equations written?
These are also given in "main_singleblock.cpp", just before the main
program starts. Non of the routines are using any property for these
equations, so they can be replaced by other ones, for instance the
nonlinear ones. One can also of course change the forcing too.