Supplementary material

to the article

Validated numerics for period-tupling and touch-and-go bifurcations of symmetric periodic orbits in reversible systems


Irmina Walawska and Daniel Wilczak

Table of content.

The purpose of this document is to provide

The C++-11 program.


The program is written in C++11. In order to compile the program one needs to

Compilation of the CAPD library.

Before you compile the CAPD library make sure that the system satisfies all requirements.

You can follow the compilation instruction or just type from the command line (under linux with gcc):

tar xvfz capd-capdDynSys-5.0.59.tar.gz
cd capd-capdDynSys-5.0.59
./configure --without-gui --with-mpfr
make -j lib

Directory structure of the program.

The main directory of the program consists of the following subdirectories

contains some configuration files, which contain numerical constants used in computation (like order of the method, tolerance, mass parameters, threshold values, etc)
contains files with precomputed data to be used by rigorous algorithms (like approximate halo orbits, approximate bifurcation points)
contains dependency files created by a compiler. They are used only to build executables
header files for all implemented general algorithms
directory used to store output of programs
contains object files
contains main programs to be invoked by the user
implementation of general algorithms

Compilation of the program.

After successful compilation of the CAPD library, enter the main directory of the program (where the makefile is located) and call

make -j CAPD=path/bin/

where path is relative or absolute path to the build directory of the CAPD library (do not miss final slash after bin). This parameter can be skipped if the capd-config file is on the system search path. Example:

make -j CAPD=/home/wilczak/capd-build/bin/

The above command builds several executables, which will be described in the next sections.

Falkner-Skan equation

The Falkner-Skan equation is given by \(f''' +ff'' - \beta((f')^2-1)=0\). The C++ program checks, that there is a smooth branch of odd periodic solutions parametrized by \(\beta\in[9/8,100\,000]\). Moreover, this branch admits period-doubling, third order touch-and-go and period-quadrupling bifurcations. All the constants (i.e. parameter range, orders of ODE solvers, approximate bifurcation points) are hard coded in a single-file prog/checkFalknerSkanBifurcations.cpp to make it as short as possible.

This program should be invoked from the command line without any parameters. One can redirect its output to a file, as shown below

./checkFalknerSkanBifurcations > log/FS.log

The program executes in a single thread within few minutes. First, it performs an adaptive subdivision of the parameter range in order to prove the existence of a smooth branch of odd periodic solutions. It logs to the standard output subinterval of parameters being processed and computed bounds on the interval Newton operator. Finally, inequalities required for the proof of three bifurcations are checked.

Output of the program: log/FS.log.

Circular Restricted Three Body Problem

The Circular Restricted Three Body Problem can be written as a system of second order equations

\[ \left\{ \begin{array}{rcl} \ddot{x} - 2\dot y &=& \frac{\partial \Omega_\mu(x,y,z)}{\partial x},\\[2pt] \ddot{y} + 2\dot x &=& \frac{\partial \Omega_\mu(x,y,z)}{\partial y},\\[2pt] \ddot{z} &=& \frac{\partial \Omega_\mu(x,y,z)}{\partial z}, \end{array} \right. \]

Numerical results regarding the above system are threefold:

Symmetry breaking bifurcations of Lyapunov orbits

Validation of symmetry-breaking bifurcation for \(L_{1,2,3}\)-halo orbits for \(\mu\in[9.5\cdot10^{-4},0.5]\) is realized by the program prog/checkSymBreakingBifurcation.cpp. It performs an adaptive subdivision of parameter range and checks the conditions of symmetry-breaking bifurcation. The program requires one parameter - a number i=1,2 or 3, which stands for the index of the family of Lyapunov orbits \(L_i\) for which computation is performed. One can redirect output of the program to a file.

./checkSymBreakingBifurcation 1 > log/symBreakingL1.log

Time of computation:
the program uses C++-11 threads to speed-up computation. The execution time varies from 7 hours for \(L_1\) to over 7 days for \(L_3\) on a computer equipped with 128 CPUs.

Data files:
The program uses precomputed database of approximate bifurcation points:

Each line consists of four numbers: mu x dy J, which stand for mass parameter, coordinates of bifurcation point and the Jacobi integral.

The program produces enormous amount of data, so we decided to store in the log file only an information about the progress of computation (ID of completed asynchronous tasks along with parameter range which has been checked). One can change this by uncommenting macro ENABLE_LOGGING in the header file include/Logger.h and recompiling the program.

Note, that the output is not deterministic due to asynchronous tasks execution. Here are sample output files:

Continuation of halo orbits in Sun-Jupiter and Earth-Moon systems

There are two programs which implement continuation of halo orbits for \(L_{1,2}\) families of halo orbits:

Both programs require a path to configuration file as a parameter. Moreover, the output can be redirected to a file. Possible usage:

./checkL1HaloContinuation conf/EarthMoonL1Branch.conf > log/EarthMoonL1Branch.log
./checkL1HaloContinuation conf/SunJupiterL1Branch.conf > log/SunJupiterL1Branch.log

./checkL2HaloContinuation conf/EarthMoonL2Branch.conf > log/EarthMoonL2Branch.log
./checkL2HaloContinuation conf/SunJupiterL2Branch.conf > log/SunJupiterL2Branch.log

The above configuration files contain hand-chosen contants and parameters of algorithms. Each configuration file contains (in particular) a path to data file with precomputed approximate selected points on branch of halo orbits. Continuation algorithm produces huge amount of data. Thus, we decided to report only progress of computation (ID of completed asynchronous tasks). We report some numbers regarding validation of symmetry breaking bifurcations and validation if computed arcs form a continuous branch.

programconfiguration filedata fileoutput
checkL1HaloContinuation conf/EarthMoonL1Branch.conf dat/Moon-L1-branch.dat log/EarthMoonL1Branch.log
checkL1HaloContinuation conf/SunJupiterL1Branch.conf dat/Jupiter-L1-branch.dat log/SunJupiterL1Branch.log
checkL2HaloContinuation conf/EarthMoonL2Branch.conf dat/Moon-L2-branch.dat log/EarthMoonL2Branch.log
checkL2HaloContinuation conf/SunJupiterL2Branch.conf dat/Jupiter-L2-branch.dat log/SunJupiterL2Branch.log

Time of computation:
The program uses C++-11 threads to speed-up computation. The execution time varies from 29 minutes for \(L_2\)-branch in the Earth-Moon system to over 29 hours for \(L_1\)-branch in the Sun-Jupiter system on a computer equipped with 128 CPUs.

Period-toupling and touch-and-go bifurcations of halo orbits.

Period-toupling and touch-and-go bifurcations in the CR3BP are validated by the program prog/checkHaloBifurcations.cpp. The program should be invoked with exactly one parameter - a path to data file.

./checkHaloBifurcations dat/filename.dat

Each data file contains the following data (in this order):

The program saves logs from computation to the file log/filename.dat_log. This file contains in particular bounds on the computed quantities as well as information if required inequalities for bifurcation are satisfied. Below we give list of links to 33 data files and corresponding output files.

Sun-Jupiter system Earth-Moon system
Input (data file) Output of the program Input (data file) Output of the program
Bifurcations of \(L_1\)-halo orbits
Bifurcations of \(L_2\)-halo orbits

Time of computation:
The program creates three asynchronous tasks (executed in three CPUs if available). Execution time of the program varies from 12 minutes to over 15 hours.

Validation of period-tupling bifurcation requires costly computation of third order derivatives of Poincare map. The configuration file contains a flag indicating, if the program should compute third order derivatives or not. This is not needed for all touch-and-go bifurcations. In the three cases of period-quadrupling bifurcations marked by star in the table above, we also set this flag to false. In these three cases, after many days of computation, the program returned bounds which were not sharp enough to conclude, that the Jacobi integral is convex (or concave) along bifurcating curve. In order to shorten time of computation we decided to skip this unsuccessful and costly step leaving validation of bifurcation incomplete (after checking that there are two intersecting curves of period-1 and period-4 points).