The PIM-triple method developed by H. Nusse and J. Yorke (see e.g., Physica D 36, 137 (1989)) allows to obtain trajectories on non-attracting sets with one-dimensional unstable manifolds. It has been widely used in the case of chaotic saddles. The PIM-simplex method is an extension of the PIM-triple method to deal with saddles that have more than one expanding direction (P. Moresco and S. Ponce Dawson, Physica D (submitted)). It combines the ideas of the PIM-triple method with an algorithm for finding local extrema of multi-variable functions (J.A Nelder and R. Mead, Computer Journal 7, 308 (1965)).
Here we describe the steps that the PIM-simplex method follows to find a non-attracting invariant set (a hyperbolic unstable periodic orbit or a chaotic saddle) with n expanding directions (an nD local unstable manifold) for a diffeomorphism of Rm into itself. In order to do this we first need to introduce some definitions.
Definitions.
A restraining region, R, is an m-dimensional subset of Rm. Given a point x in R we define its escape time, T(x), as the minimum integer N>0 such that fN(x) is not in R. Given a point x in R with finite escape time, we define the escape distance, d(x), as the sum over i from i=0 to i=T(x) of |fi(x) - fi-1(x) | . An n-dimensional simplex is the geometrical figure consisting of n+1 points (or vertices), P1, P2,...,Pn+1, and all their interconnecting line segments, polygonal faces, etc. In two dimensions, it is a triangle, in three dimensions, a tetrahedron (not necessarily regular), etc. We will only consider non-degenerate simplices, i.e., those that enclose a finite inner n-dimensional volume (W.H. Press et al, "Numerical Recipes", Cambridge University Press, 1992, p. 408).
The method.
1. A restraining region, R, is defined. It need not be singly connected. The invariant set that one is looking for needs to be entirely contained in R.
2. n+1 points, P1, P2,...,Pn+1, defining an n-dimensional simplex, S0, are chosen in a way such that all segments connecting P1 to the other points are orthogonal among themselves and of size (much) bigger than a tolerance d1.
3. The Simplex method of Nelder and Mead is applied to S0 in order to find a maximum of the escape time. If during this procedure a situation is reached such that all points defining the simplex have the same escape time and the length of the simplex sides is bigger than a tolerance d0 < d1, the escape distance is evaluated, instead, to make transformations on the simplex. In this case the transformation is chosen so as to decrease the escape distance.
4. Step 3 is repeated until a simplex, S0', is found such that the length of all its sides is less than d0<<d1 (d0 is usually taken very close to machine precision).
5. A new simplex, S0'', is chosen, entirely contained
in S0, defined by the following n+1 points:
the point in S0' with the largest escape time, Pmax,
and n other points, a distance d1 apart from Pmax,
so that the segments they form with Pmax are at right angles
among themselves. If d1 is small enough, S0''
provides an approximation for a point, x0, on the stable manifold
of the invariant set we are looking for.
6. The map is applied to the n+1 points that define S0'' and then the steps 3--4 are repeated. In this way we get an approximation for f(x0).
7. By repeated applications of steps 3--6 we get an approximation for
a trajectory that starts close to the stable manifold of the saddle that
eventually approaches the saddle.
Examples and credits.
This method was developed by Pablo
Moresco and Silvina Ponce Dawson. We benefited from discussions and
suggestions from Gabriel Mindlin and Jim Yorke. Pablo did all the
coding which he is now making public for anybody to use (and modify). All
the necessary files are tared up in a single file which we called pims.tar.
Once you retrieve it you need to untar it (if you
save it as pims.tar, type: tar -xvf pims.tar). When you do
that, a directory called pims-code containing
all these files will be created. The code is written in C and runs fine
on a Linux platform. We have not checked other platforms. The version we
are making public contains three examples:
the unstable fixed point of a very simple map, a chaotic saddle of the
double-rotor map (see e.g., C. Grebogi, E. Kostelich, E. Ott, and
J. A. Yorke, Physica D 25, 347 (1987)) and a chaotic saddle of a
three-dimensional horseshoe map. You choose the example in the Makefile
(see the README file for further explanations). In all cases we use d0=10-12
and d1=10-8 and the parameters suggested
by Nelder and Mead in their original paper. The example for the double-rotor
map is a little bit tricky and can only run for 600 iterates or so. If
you want to go beyond that, you need to divide the restraining region in
a lot of singly connected pieces (not only two, as it is now in the code).
Those interested in looking at that should contact Pablo (email:
mad16@cc.keele.ac.uk).