A Hamiltonian system of three degrees of freedom with eight channels of escape: The Great Escape
Abstract
In this work, we try to shed some light to the nature of orbits in a threedimensional potential of a perturbed harmonic oscillator with eight possible channels of escape, which was chosen as an interesting example of open threedimensional Hamiltonian systems. In particular, we conduct a thorough numerical investigation distinguishing between regular and chaotic orbits as well as between trapped and escaping orbits, considering unbounded motion for several values of the energy. In an attempt to discriminate safely and with certainty between ordered and chaotic motion, we use the Smaller ALingment Index (SALI) detector, computed by integrating numerically the basic equations of motion as well as the variational equations. Of particular interest, is to locate the basins of escape towards the different escape channels and connect them with the corresponding escape periods of the orbits. We split our study into three different cases depending on the initial value of the coordinate which was used for launching the test particles. We found, that when the orbits are started very close to the primary plane the respective grids exhibit a high degree of fractalization, while on the other hand for orbits with relatively high values of several wellformed basins of escape emerge thus, reducing significantly the fractalization of the grids. It was also observed, that for values of energy very close to the escape energy the escape times of orbits are large, while for energy levels much higher than the escape energy the vast majority of orbits escape extremely fast or even immediately to infinity. We hope our outcomes to be useful for a further understanding of the escape process in open 3D Hamiltonian systems.
Keywords:
Hamiltonian systems; harmonic oscillators; numerical simulations; escapes; fractals∎
1 Introduction
The issue of escapes in Hamiltonian systems is directly related to the problem of chaotic scattering which has been an active field of research over the last decades and it still remains open (e.g., [BTS96, – BGOB88, , CPR75, , C90, , CK92, , E88, , JS88, , ML02, – PH86, , SASL06, – SS10, ]). We know, that particular types of Hamiltonian systems have a finite energy of escape and for energy levels beyond the escape energy, these systems allow particles to escape to infinity. The literature is replete with studies of such “open” Hamiltonian systems (e.g., [BBS09, , CKK93, , KSCD99, , SCK95, – SKCD96, ]).
Usually, the infinity acts as an attractor for an escape particle, which may escape through different channels (exits) on the equipotential curve or on the equipotential surface depending whether the dynamical system is two or threedimensional, respectively. Therefore, it is quite possible to obtain basins of escape, similar to basins of attraction in dissipative systems or even the NewtonRaphson fractal structures. Basins of escape have been studied in several papers (e.g., [BGOB88, , C02, , KY91, , PCOG96, ]). The reader can find more details regarding basins of escape in [C02, ].
The wellknown HénonHeiles system [HH64, ] is undoubtedly a paradigmatic model for timeindependent Hamiltonian systems of two degrees of freedom. A huge load of research on escape from this system has been conducted over the years (e.g., [AJ03, – AVS03, , BBS08, , BSBS12, , dML99, , S07, ]). Here, we would like to point out that all the above references on escapes in the HénonHeiles system are exemplary rather than exhaustive, taking into account that a vast quantity of related literature exists.
A simple dynamical system of two coupled harmonic oscillators for various values of the energy above the escape energy has been investigated in [CE04, ], where it was found that stable periodic orbits are surrounded by stability islands that never escape. A further numerical analysis of the same dynamical system in [CHLG12, ], revealed that as the energy increases beyond the escape energy, the majority of chaotic orbits escape either directly, or after a small or large number of intersections with the axis. In the same vein, the effects of different types of perturbations on both the topology and the escaping dynamics in the HénonHeiles system was studied in [BSBS12, ], where basins of escape were found to exits in the physical as well as in the phase space.
Of particular interest is the issue of escaping orbits in galactic dynamics. In a recent article [Z12b, ], we explored the nature of the orbits of stars in a galactictype potential, which can be considered to describe local motion in the meridional plane near the central parts of an axially symmetric galaxy. It was observed, that apart from the trapped orbits there are two types of escaping orbits, those which escape fast and those which need to spend vast time intervals inside the limiting curve before they find the exit and eventually escape. Furthermore, the chaotic dynamics within a star cluster embedded in the tidal field of a galaxy was examined in [EJSP08, ]. In particular, by scanning thoroughly the phase plane and obtaining the basins of escape with the respective escape times it was revealed, that the higher escape times correspond to initial conditions of orbits near the fractal basin boundaries.
We should like to point out, that the vast majority of the existed literature deals with escaping orbits in Hamiltonian systems of two degrees of freedom, while only a handful of papers is devoted to escapes in dynamical systems with three degrees of freedom. Thus, we decided to make a new contribution by exploring the nature of orbits and the escape process in a simple threedimensional (3D) potential with eight possible channels of escape. The aim of this work, is twofold: (i) to distinguish between ordered/chaotic and trapped/escaping orbits and (ii) to locate the basins of escape leading to different escape channels and try to connect them with the corresponding escape times of the orbits.
Over the last half century, dynamical systems made up of perturbed harmonic oscillators have been extensively used in order to describe local motion (i.e., near an equilibrium point) (e.g., [AEFR06, , C93, , CK98, , CK99, , FLP98a, – HH64, , SI79, , Z12c, ]). In an attempt to reveal and understand the nature of orbits in these systems, scientists have used either numerical (e.g., [KV08, , ZC12, ]) or analytical methods (e.g., [CB82, , D91, , DE91, , E00, , ED99, ]). Furthermore, potentials made up of harmonic oscillators are frequently used in Astronomy, as a first step for distinguishing between ordered and chaotic local motion in galaxies, since it is widely accepted that the motion of stars near the central region of a galaxy can be approximated by harmonic oscillations. Therefore, the abovementioned facts justify our choice of using a potential of a perturbed harmonic oscillator in our quest for escaping orbits.
The structure of the present article is as follows: in Section 2 we describe the properties of the potential we chose for our investigation of escaping orbits. The computational methods used in order to determine the nature (ordered/chaotic and trapped/escaping) of orbits are described in Section 3. In the following Section, we conduct a thorough analysis of the orbits presenting in detail all the numerical results of our computations. Our article ends with Section 5, where the discussion and the conclusions of this research are presented.
2 Properties of the model potential
The general form of a threedimensional perturbed harmonic oscillators is
(1) 
where , and are the unperturbed frequencies of oscillations along the , and axes respectively, is the perturbation parameter, while is the function containing the perturbing terms. This is called a threedimensional perturbed elliptic oscillator.
In the present paper, we shall use a threedimensional perturbed harmonic oscillator at the 1:1:1 resonance, that is when , in order to investigate the escape properties of orbits. The corresponding potential is
(2) 
where in our case the perturbation is given by
(3) 
being the common frequency of oscillations along the three axes. Without the loss of generality, we may set and for more convenient numerical computations. The same potential was also used in [CZ12, ] however, the study of the nature of orbits was restricted only to bounded motion.
The Hamiltonian to potential (2) reads
(4) 
where , and are the momenta per unit mass conjugate to , and respectively, while is the numerical value of the Hamiltonian, which is conserved.
The basic equations of motion for a test particle with a unit mass are
(5) 
where, as usual, the dot indicates derivative with respect to the time. Furthermore, the variational equations governing the evolution of a deviation vector are
Potential (2) has a finite energy of escape which is given by
(7) 
where . Eq. (7) is a general formula, while in our case where the escape energy takes the value . For values of energy smaller than the energy of escape, the isopotential surface is closed. We see in Fig. 1a where , that the threedimensional isopotential surface is closed and has a 3D starlike shape. In Fig. 1b where we observe the presence of some very sharp edges indicating that the surface is about to open. Indeed, when it can be seen in Fig. 1c, that the surface opens and eight holes appear. These holes are in fact escaping channels through which a test particle is free to escape to infinity only when the motion is unbounded . The escaping channels are symmetrical to the three primary planes , due to the overall symmetry of the surface. In our investigation, we will perform a thorough statistical analysis regarding the escaping channels therefore, in Fig. 2 we number each channel. Even though the escape channels are symmetrical, we can distinguish between them from their particular position on the surface and the corresponding sings of the coordinates . In particular we have: , , , , , , and , where indicates positive coordinate, while stands for negative coordinates.
An issue of paramount importance is the determination of the position as well as the time at which an orbit escapes. In twodimensional systems, the highly unstable Lyapunov periodic orbits control the properties of escaping orbits as they bridge the openings of the zero velocity curve. In particular, any orbit that crosses a Lyapunov orbit with velocity pointing outwards escapes from the system (see e.g., [C90, ]). As far as we know, in threedimensional systems this method cannot be applied so, we have to find another, possibly new, way of determining the escapes of orbits. Taking into account the particular geometry of the isopotential surface, we could define an appropriate cutoff surface which embraces the isopotential and abuts on the eight escape channels (see Fig. 3). After conducting extensive numerical simulations we found, that the best cutoff surface is the following
(8) 
Thus, a threedimensional orbit is considered to escape one when . In Fig. 4 we present a characteristic example of an orbit which escapes from channel 1.
In our investigation, we shall deal only with unbounded motion of test particles for values of energy in the set . Our numerical calculations indicate, that when the value of the energy is large enough the openings are so wide, that as we can see in Fig. 1d where the isopotential surface is dismembered into six distinct components and therefore, it is no longer possible to define escapes through particular channels.
3 Computational methods
Knowing the ordered or chaotic nature of orbits (trapped and escaping) is an issue of great interest. Over the years, several chaos indicators have been developed in order to determine the character of orbits. In our case, we chose to use the Smaller ALingment index (SALI) method. The SALI [S01, , SABV04, ] has been proved a very fast, reliable and effective tool, which is defined as
(9) 
where and are the alignments indices, while and , are two deviations vectors which initially point in two random directions. For distinguishing between ordered and chaotic motion, all we have to do is to compute the SALI for a relatively short time interval of numerical integration . In particular, we track simultaneously the timeevolution of the main orbit itself as well as the two deviation vector and in order to compute the SALI. The variational equations (LABEL:variac), as usual, are used for the evolution and computation of the deviation vectors.
The timeevolution of SALI strongly depends on the nature of the computed orbit since when the orbit is regular the SALI exhibits small fluctuations around non zero values, while on the other hand, in the case of chaotic orbits the SALI after a small transient period it tends exponentially to zero approaching the limit of the accuracy of the computer . Therefore, the particular timeevolution of the SALI allow us to distinguish fast and safely between regular and chaotic motion. The timeevolution of a regular (R) and a chaotic (C) orbit for a time period of time units is presented in Fig. 5a. We observe, that both regular and chaotic orbits exhibit the expected behavior. Nevertheless, we have to define a specific numerical threshold value for determining the transition from regularity to chaos. After conducting extensive numerical experiments, integrating many sets of orbits, we conclude that a safe threshold value for the SALI is the value . The horizontal, blue, dashed line in Fig. 5a corresponds to that threshold value which separates regular from chaotic motion. In order to decide whether an orbit is regular or chaotic, one may use the usual method according to which we check after a certain and predefined time interval of numerical integration, if the value of SALI has become less than the established threshold value. Therefore, if SALI the orbit is chaotic, while if SALI the orbit is regular. In Fig. 5b we present a dense grid of initial conditions , where the initial conditions are colored according to the value of SALI. We clearly observe several regions of regularity indicated by light reddish colors as well as a unified chaotic domain (blue/purple dots). Therefore, the distinction between regular and chaotic motion is clear and beyond any doubt when using the SALI method.
A simple qualitative way for distinguishing between regular and chaotic motion in a Hamiltonian system is by plotting the successive intersections of the orbits using a Poincaré Surface of Section (PSS) [HH64, ]. This method has been extensively applied to 2D models, as in these systems the PSS is a twodimensional plane. In 3D systems, however, the PSS is fourdimensional and thus the behavior of the orbits cannot be easily visualized. One way to overcome this issue is to project the PSS to phase spaces with lower dimensions, following the method used in [Z12a, , Z12c, , ZC13, ]. Let us start with initial conditions on a 4D grid of the PSS. In this way, we are able to identify again regions of order and chaos, which may be visualized, if we restrict our investigation to a subspace of the whole 6D phase space. We consider orbits with initial conditions , , while the initial value of is always obtained from the energy integral (4). In particular, we define a value of , which is kept constant and then we calculate the SALI of the 3D orbits with initial conditions , . Thus, we are able to construct again a 2D plot depicting the plane but with an additional value of , since we deal with 3D motion. All the initial conditions of the 3D orbits lie inside the limiting curve defined by
(10) 
For the study of our models, we need to define the sample of orbits whose properties (chaos or regularity, escaping or trapped) we will identify. The best method for this purpose, would have been to choose the sets of initial conditions of the orbits from a distribution function of the model. This, however, is not available so, we define, for each set of values of the energy, a dense grid of initial conditions regularly distributed in the area allowed by the value of the energy, following the abovementioned method. In each grid the step separation of the initial conditions along the and axis was controlled in such a way that always there are at least 25000 orbits. For each initial condition, we integrated the equations of motion (5) as well as the variational equations (LABEL:variac) using a double precision BulirschStoer FORTRAN algorithm (e.g., [PTVF92, ]) with a small time step of order of , which is sufficient enough for the desired accuracy of our computations (i.e. our results practically do not change by halving the time step). In all cases, the energy integral (Eq. (4)) was conserved better than one part in , although for most orbits it was better than one part in .
In our computations, we set time units as a maximum time of numerical integration. The vast majority of orbits (regular and chaotic) however, need considerable less time to find one of the exits in the limiting surface and eventually escape from the system (obviously, the numerical integration is effectively ended when an orbit passes through one of the escape channels and intersects the cutoff surface). Nevertheless, we decided to use such a vast integration time just to be sure that all orbits have enough time in order to escape. Remember, that there are the so called “sticky orbits” which behave as regular ones during long periods of time. A characteristic example of a sticky orbit (S) in our dynamical model can be seen in Fig. 5a, where we observe that the chaotic character of the particular sticky orbit is revealed only after a considerable long integration time of about 10000 time units. Here we should clarify, that orbits which do not escape after a numerical integration of time units are considered as nonescaping or trapped.
4 Numerical results
Our main objective is to determine which orbits escape and which remain trapped, distinguishing simultaneously between regular and chaotic motion. Moreover, two additional properties of the orbits will be examined: (i) the channels through which the particles escape and (ii) the timescale of the escapes (we shall also use the term escape period). In the present paper, we explore these aspects for various values of the energy , as well as for the initial value of the coordinate. In particular, three different cases are considered: (a) orbits starting very close to the plane (low value of ), (b) orbits with a mediocre value of and (c) orbits launched at relatively large distance from the plane (high value of ).
4.1 Case I: Orbits starting with a low value of
We start our exploration, considering orbits staring very close to the primary plane with initial condition . Fig. 6 shows the structure of the plane for several values of the energy, where different colors indicate the four types of orbits. Specifically, black color corresponds to regular trapped orbits, blue color to chaotic trapped orbits, green color to regular escaping orbits, while the initial conditions of chaotic escaping orbits are marked with red color. The outermost black solid line is the limiting curve defined by Eq. (10). We observe, that for all studied energy levels the vast majority of the initial conditions correspond to either regular or chaotic escaping orbits, while trapped orbits are mainly confined to the center of the grid where small stability regions should exist. Here we should notice, that generally there is no indication of welldefined regions corresponding to a particular type of orbits. In fact, the initial conditions of the four types of orbits seem to be randomly scattered across the plane. It is also seen, that for low values of the energy (very close to the energy of escape) almost the entire grid is covered by chaotic escaping orbits. As the energy increases however, regular escaping orbits are gaining ground and finally they take over the grid. All these can be seen better in Fig. 7 where we present the evolution of the percentages of regular and chaotic escaping orbits as a function of the energy . Indeed, the rates evolve similarly but following a complete different direction; chaotic escaping orbits decrease, while regular escaping orbits increase. It is worth noticing, that about they share the whole area of the grid, while at the highest studied value of energy regular escaping orbits are four time more the chaotic escaping orbits^{1}^{1}1Generally, any dynamical method requires a sufficient time interval of numerical integration in order to distinguish safely between ordered and chaotic motion. Therefore, if the escape rate of orbits is very low or even worse if the orbits escape directly from the system then, any chaos indicator (the SALI in our case) will fail to work properly due to insufficient integration time. Nevertheless, we decided to apply the SALI method regardless of the escape rate of orbits..
Fig. 8 depicts in another point of view the structure of the plane. In this case, each initial condition is colored according to the escape channel through which the particular orbit escapes. The black regions denote initial conditions where the test particles do not escape. Looking those grids we can draw two conclusions: (i) there is a highly sensitive dependence of the escape process on the initial conditions, that is, a slight change in the initial conditions makes the test particle escape through another channel, which is is a classical indication of chaos; (ii) All computed grids exhibit a very rich fractal structure regarding the escaping channels. This is true, because there is no indication whatsoever of basins of escape. By the term basin of escape, we refer to a set of initial conditions that corresponds to a certain escape channel. Inspecting carefully the grids shown in Fig. 8, one may identify very small and local patterns leading to groups of neighbour initial conditions with the same color, however, these groups are extremely small and therefore cannot be recognized as basins of escape. The fractalization of the grid indicates that all channels of escape should be equiprobable. This is confirmed in Fig. 9 where we present a histogram of the percentages of escaping orbits per channel for all the tested values of the energy. It is seen, that in general terms (in some cases, the rates of channels 7 and 8 are slightly elevated) all channels have the same percentage, roughly around 12.5% (indicated by the horizontal black dashed line), which yields to equal probability of 1/8. In Fig. 10 we present the isopotential surface for , where we have marked with red dots the escape positions^{2}^{2}2The term “escape positions” refers to the points of the 3D space at which the orbits intersect the cutoff surface with velocity pointing outwards and escape. of the orbits. Taking into account that the corresponding grid contains a very rich sample of orbits (about 25000) and almost all of them escape choosing randomly and therefore equiprobably a channel, we see that all eight exits of the surface have been patched equally well thus verifying that in this case the escape channels are indeed equiprobable.
The following Fig. 11 shows how the escape times of orbits are distributed on the grid. Light reddish colors correspond to fast escaping orbits, dark blue/purpe colors indicate large escape periods, while black color denote trapped orbits. We observe, that when , that is a value of energy very close to the escape energy, the escape periods of the majority of orbits are huge corresponding to tens of thousands of time units. This however, is anticipated because in this case the width of the escape channels is very small and therefore, the orbits should spend much time inside the isopotential surface until they find one of the openings. As the value of the energy increases, the escape channels become more and more wide leading to faster escaping orbits, which means that the escape period decreases rapidly. To prove this point, we calculated for each value of the energy the respective average value of the escape period . Our results are given in Fig. 12, where we see that for the escape period of orbits decreases rapidly, while for larger energies it follows an almost linear decrease. Here we should note, that for large values of the energy weak color patterns of elongated bands spiralling around the center of the grids appear. These patterns correspond to extremely fast escaping orbits with escape periods less than 10 time units.
4.2 Case II: Orbits starting with a mediocre value of
Our investigation continues considering in this case orbits starting at a moderate distance from the primary plane having initial condition . We shall follow the same numerical approach as in Case I. The structure of the plane using different colors in order to distinguish between the four types of orbits (regular trapped; chaotic trapped; regular escaping and chaotic escaping) is presented in Fig. 13. We observe similar results to that discussed previously in Fig. 6. It is worth noticing, that the increase of the initial value of the coordinate of orbits, has affected the geometry of the grid since it has become more oval. Our numerical calculations reveal, that in most examined cases all the integrated orbits of the grids escape sooner or later, while only the periodic orbits do not escape and only when they are stable^{3}^{3}3The stability of a periodic orbit is determined by computing the stability index [Z13, ]. A periodic orbits is stable if only the stability index (S.I.) is between 2 and +2.. We also see, that for low values of the energy chaotic escaping orbits swarm the grid, while as the energy increases regular escaping orbits become the most populated type of orbits. Fig. 14 shows how the energy influence the percentages of regular and chaotic escaping orbits. The pattern of the evolution is quite similar to that of Fig. 7 however, there are two noticeable differences: (i) the point where regular and chaotic escaping orbits share the grid comes earlier at about and (ii) at the highest studied energy level regular escaping orbits are nine time more the chaotic escaping orbits.
The plane is reconstructed in Fig. 15, where each initial condition is colored according to the escape channel through which the particular orbit escapes. It is seen, that for low values of the energy the grids exhibit a large degree of fractalization regarding the escape channels. However, when small groups of neighbor initial conditions leading to the same exit start to emerge. As the energy increases, these structures expand in size thus forming welldefined basins of escape, while the degree of fractalization of the grid is consequently reduced. The basins of escape are either broad regions, or elongated bands spiralling around the center of the grid. The various basins of escape are very intricately interwoven. The formation of several basins of escape implies that the escape channels should not be equiprobable any longer. This assumption is justified in Fig. 16, where a histogram of the percentages of escaping orbits per channel for all the tested values of the energy is given. One can observe, that when the test particles seem not to choose randomly an exit but they prefer specific channels through which they escape. In general terms, channels 1, 2, 7 and 8 appear to be the most favorite ones and are more likely to be chosen, thus having a probability higher than the average. Exits 3 and 6 on the other hand, are the least favorite having the smallest probability to be chosen.
The distribution of the escape times of orbits on the grid when is shown in Fig. 17. We observe several similarities with Fig. 11 which discussed in the previous case where being the most noticeable the large escape periods of orbits for values of the energy very close to the escape energy and the rapid decrease of the escape times with increasing energy. This behavior can be comprehended better by looking at Fig. 18 where the evolution of the average value of the escape period of orbits as a function of the energy is given. We see, that the average escape time displays a rapid reduction when , while for larger energy values it decreases almost linearly. At this point, we should notice that the basins of escape for , can be distinguished in the grids of Fig. 17, being the regions with intense red color indicating extremely fast escaping orbits. Our numerical calculations indicate, that orbits with initial conditions inside the basins have significantly small escape periods of less than 5 time units.
4.3 Case III: Orbits starting with a high value of
The last case under investigation involves orbits launched at a relatively high distance from the plane having initial condition . Again, the numerical approach remains the same as in the previously studied cases. Fig. 19 shows the structure of the plane distinguishing between ordered/chaotic and trapped/escaping orbits using different colors. We see that things are quite similar to those of Figs. 6 and 13 however, the shape of the grid, or in other words the liming curve, has become much more oval and elongated along the axis, due to the further increase of the initial value of the coordinate. As in the previous cases, only the central periodic orbit at remains trapped when it is stable, while all the other orbits escape from the system. Once more it is evident, that at low values of the energy the majority of escaping orbits are chaotic, while at high values of the energy regular escaping orbits dominate the grid. The evolution of the percentages of regular and chaotic escaping orbits as a function of the energy is given in Fig. 20. It can be seen, that at about regular and chaotic escaping orbit share the entire grid, while at high energy levels , chaotic escaping orbits are depleted and more than 95% of the plane is covered by fast regular escaping orbits.
In Fig. 21 we present another point of view of the grids shown in Fig. 19. Now, the color code of the initial conditions of orbits depends on the escape channels. One may observe, that apart from the gird in which the structure of the plane displays high fractalization, in all other studied cases several basins of escape are present. At low values of the energy, we observe that the basins of escape appear only as weak elongated bands spiralling around the center, while as we proceed to grids corresponding to higher energy levels broad escape regions emerge which eventually take over the majority of the grid at high energies . Therefore, we may conclude that the increase of the energy favors the growth and the development of basins of escape. Fig. 22 shows the percentages of escaping orbits per channel for all the tested values of the energy. Evidently, the escape channels can be considered equiprobable only at values of energy very close to the escape energy . In particular, for channels 3 and 4 seem to be the most favorable, while this behavior is reversed for larger values of the energy. Indeed, we observe that when channels 1 and 2 are clearly the most selected ones, while channels 3 and 4 are the least chosen.
Finally, Fig. 23 shows the distribution of the escape times of orbits on the grids when . The results are very similar to those presented in Figs. 11 and 17. We see, that the general rule according to which the escape times decrease with increasing energy is still valid. The evolution of the average value of the escape period of orbits as a function of the energy given in Fig. 24 is also in the same vein. As we pointed out in the previous case, the several basins of escape are visible in the grids of Fig. 23, being the regions with intense red color indicating orbits with extremely low escape periods. According to our numerical computations, orbits with initial conditions inside the basins escape almost immediately from the system.
4.4 General overview
Before closing this section, we would like to summarize and compare the results regarding the basins of escape. In Fig. 25 we show the evolution of the percentage of the total area on the plane corresponding to basins of escape, either broad regions or elongated spiral bands, as a function of the energy for the three initial values of the coordinate. It is seen, that for low values of (see Case I in subsection 4.1, where ), only at high energy levels we have weak indications of basins of escape corresponding to a tiny fraction of the grid, less than 5%, while the vast majority of the grid displays a highly fractal structure. When on the other hand, basins of escape appear quite early for and for larger values of the energy, they start to grow rapidly occupying about 65% of the grid at the highest studied energy level. The growth of the extent of escape basins is considerably faster in the case where orbits are launched with high values of . We observe, that for the rate of the basins of escape increases almost linearly obtaining high values very fast and covering about 90% of the entire grid when . Therefore, one may conclude that for low values of the degree of fractalization is very high, while as we proceed to larger values of this behavior is restricted or confined and wellformed basins of escape take over the majority of the planes.
One may reasonably assume, that the basins of escape appear only at the plane and if we choose another subspace of the entire 6D space they might not be present. In order to give an answer to this, we decided to integrate orbits in three more planes and therefore check if the basins of escape are still present or not. Our results are presented in Fig. 26(ad) where the model is examined. In particular, Fig. 26a shows a grid of initial condition in the physical plane when , , Fig. 26b depicts the plane discussed previously in Fig. 21, Fig. 26c shows a grid of initial condition in the plane when , , , while in Fig. 26d the grid of initial condition in the plane when and , while in all cases is obtained from the energy integral (4). It is evident, that the escape basins are indeed present in all studied cases thus implying that their formation is an intrinsic property of the dynamical system regardless of the particular chosen twodimensional plane of the 6D phase space.
Our numerical analysis reveals, that orbits with initial conditions inside the basins escape have either extremely low escape periods or they escape immediately from the system. At this point, we would like to emphasize this point a little bit. For this purpose, we constructed a threedimensional plot shown in Fig. 27, which combines results regarding both the escape channels and escape rates of the model. It is seen, that obits belonging to escape basins possess low escape rates, while on the other hand, orbits inside the fractal region have relatively larger escape rates. With a much closer look we can see that orbits inside broad escape regions escape almost immediately, while orbits inside elongated bands spiralling around the center escape a little bit later.
So far in our investigation, we chose discrete energy levels and we explored the escape process of orbits for three different values of the coordinate. It would be of particular interest, to study how the basins of escape evolve within a dense spectrum of values of the energy when . Therefore, we integrated orbits with initial conditions on the axis with , , while the value of is obtained from the energy integral. In Fig. 28a we distinguish the initial conditions of orbits with respect to the escape channel, while Fig. 28b shows the plane with the corresponding escape periods. The outermost black solid line denotes the Zero Velocity Curve (ZVC). We observe, that as we proceed to higher energy levels the presence of the escape basins becomes more and more prominent and at the same time the escape period reduces, while on the other hand the central region of the plane remains fractal all the way up. We would also like to point out, that both plots are qualitatively symmetrical with respect to the axis.
Our last step will be to try to consolidate a general view of the escape process for different values of the coordinate. To obtain this, we defined a dense, threedimensional grid of about initial conditions with inside the isopotential surface when . As usual, the initial value of was derived from the energy integral. In Fig. 29 we visualize the results of the numerical integration showing with the different color the escape channels. One may observe, once more, several basins of escape embedded in the surface of the grid. In this case, the grid is a threedimensional solid therefore, only its surface is visible. However, we can penetrate inside the solid by using a tomographicstyle approach. According to this method, we can again plot twodimensional slices of the solid grid by defying specific levels to each primary plane. Fig. 30 shows the evolution of the structure of the primary , and planes when . We see, that the structure evolves rapidly and non uniformly (constant interplay between escape basins and fractal structure) thus implying that the escape process in this Hamiltonian system with three degrees of freedom (3 d.o.f) is, by all means, a very complex procedure.
5 Discussion and conclusions
In the present paper, we tried to shed some light to the escape process in threedimensional open Hamiltonian systems. We chose a simple potential of a perturbed harmonic oscillator which has a finite energy of escape. Beyond the escape energy the threedimensional isopotential surfaces open and eight symmetric channels of escape appear. Our main objective was to reveal the nature of orbits and try to understand the escape process. In particular, we managed to distinguish between ordered/chaotic and trapped/escaping orbits and we also located the basins of escape leading to different escape channels, finding correlations with the corresponding escape times of the orbits. Our extensive and thorough numerical investigation strongly suggests, that overall escape process is very dependent on the value of the total orbital energy, as well as on the initial value of the coordinate.
Since a distribution function of the model was not available so as to use it for extracting the different samples of orbits, we had to follow an alternative path. We decided to restrict our exploration to a subspace (a 4D grid) of the whole 6D phase space. Specifically, we considered orbits with initial conditions , , while the initial value of is always obtained from the energy integral. Then, we defined a value of , which was kept constant. Thus, we were able to construct again a 2D plot depicting the plane (we may call it a local phase plane) but with an additional value of , since we deal with 3D motion. Following this approach, we were able to identify again regions of order and chaos. We defined for several values of the energy , a dense grid of initial conditions regularly distributed in the area allowed by the value of the energy. In each grid the step separation of the initial conditions along the and axis was controlled in such a way that always there are at least 25000 orbits. For the numerical integration of the orbits in each grid, we needed about between 2 minutes and 5.1 hours of CPU time on a Pentium DualCore 2.2 GHz PC, depending on the escape rates of orbits in each case.
In this work, we revealed the influence of the orbital energy , as well as the initial value of the coordinate of threedimensional orbits on the escape process. The main numerical results of our research can be summarized as follows:

It was found, that more than 95% of the computed orbits escape sooner or later, while only a tiny fraction of orbits corresponding mainly to stable periodic orbits remain trapped. For values of energy very close to the escape energy, the majority of escaping orbits are chaotic, while the percentage of chaotic escaping orbits reduces rapidly with increasing energy and for high energy levels the ordered escaping orbits is the most populated type of orbits. This behavior is more or less the same, regardless the particular initial value of the coordinate.

A strong correlation between the extent of the basins of escape and was found to exists. Our numerical calculations indicate that for low values of the structure of the phase plane exhibits a large degree of fractalization and therefore orbits escape choosing randomly escape channels. As the value of increases however, the grid becomes less and less fractal and several basins of escape, either broad regions or elongated bands spiralling around the center emerge. The extent of these basins of escape is more prominent at high values of the energy, where they occupy about 90% of the entire grid.

It was also observed, that the escape times of orbits are directly linked to the basins of escape. Specifically, for low values of energy close to , the escape periods are huge corresponding to tens of thousands of time units. This was justified, by taking into account that for these energy levels the width of the escape channels is very small, thus forcing the orbits to consume large time intervals inside the isopotential surface until they locate an exit and eventually escape to infinity. For higher energy levels on the other hand, several basins of escape are present so, the orbits escape very quickly within less than 10 time units, or even immediately from the system.
The outcomes and the conclusions of the present research are considered, as an initial effort and also as a promising step in the task of understanding the escape mechanism of orbits in open threedimensional Hamiltonian systems. Taking into account that our results are encouraging, it is in our future plans to expand our investigation in other more complicated potentials, focusing our interest in reveling the escape process of stars in galactic systems such as star clusters, or binary stellar systems.
Acknowledgments
The author would like to thank the two anonymous referees for the careful reading of the manuscript and for all the apt suggestions and comments which allowed us to improve both the quality and the clarity of the paper.
References
 (1) Aguirre, J., Sanjuán, M.A.F.: Limit of small exits in open Hamiltonian systems. Phys. Rev. E 67, 0562011–7 (2003)
 (2) Aguirre, J., Vallego, J.C., Sanjuán, M.A.F.: Wada basins and chaotic invariant sets in the HénonHeiles system. Phys. Rev E 64, 0662081–11 (2001)
 (3) Aguirre, J., Vallego, J.C., Sanjuán, M.A.F.: Wada basins and unpredictability in Hamiltonian and dissipative systems. Int. J. Mod. Phys. B 17, 41714175 (2003)
 (4) Arribas, M., Elipe, A., Floria, L., Riaguas, A.: Oscillators in resonance . Chaos Solitons Fractals 27, 12201228 (2006)
 (5) Barrio, R., Blesa, F., Serrano, S.: Fractal structures in the HénonHeiles Hamiltonian. Europhys. Lett. 82, 100031–6 (2008)
 (6) Barrio, R., Blesa, F., Serrano, S.: Bifurcations and safe regions in open Hamiltonians. New J. Phys. 11, 0530041–12 (2009)
 (7) Benet, L., Trautman, D., Seligman, T.: Chaotic scattering in the Restricted ThreeBody Problem. I. The Copenhagen Problem. Celest. Mech. Dyn. Astron. 66, 203228 (1996)
 (8) Benet, L., Seligman, T., Trautman, D.: Chaotic scattering in the Restricted ThreeBody Problem II. Small mass parameters. Celest. Mech. Dyn. Astron. 71, 167189 (1998)
 (9) Bleher, S., Grebogi, C., Ott, E., Brown, R.: Fractal boundaries for exit in Hamiltonian dynamics. Phys. Rev. A 38, 930938 (1988)
 (10) Blesa, F., Seoane, J.M., Barrio, R., Sanjuán, M.A.F.: To escape or not to escape, that is the question  Perturbing the HénonHeiles Hamiltonian. Int. J. Bifurc. Chaos 22, 12300101–9 (2012)
 (11) Caranicolas, N.D.: The 1:1 resonance in galactic type Hamiltonian systems. Astron. Astrophys. 267, 388391 (1993)
 (12) Caranicolas, N.D., Barbanis, B.: Periodic orbits in nearly axisymmetric stellar systems. Astron. Astrophys. 114, 360366 (1982)
 (13) Caranicolas, N.D., Karanis, G.I.: Chaos in barred galaxy models. Astrophys. Space Sci. 259, 4556 (1998)
 (14) Caranicolas, N.D., Karanis, G.I.: Motion in a potential creating a weak bar structure. Astron. Astrophys. 342, 389394 (1999)
 (15) Caranicolas, N.D., Zotos, E.E.: Investigating the nature of motion in 3D perturbed elliptic oscillators displaying exact periodic orbits. Nonlin. Dyn. 69, 17951805
 (16) Churchill, R., Pecelli, G., Rod, D.: Isolated unstable periodic orbits. J. Differ. Equ. 17, 329348 (1975)
 (17) Contopoulos, G.: Asymptotic curves and escapes in Hamiltonian systems. Astron. Astrophys. 231, 4155 (1990)
 (18) Contopoulos, G.: Order and Chaos in Dynamical Astronomy. Springer, Berlin (2002)
 (19) Contopoulos, G., Efstathiou, K.: Escapes and recurrence in a simple Hamiltonian system. Celest. Mech. Dyn. Astron. 88, 163183 (2004)
 (20) Contopoulos, G., Kaufmann, D.: Types of escapes in a simple Hamiltonian system. Astron. Astrophys. 253, 379388 (1992)
 (21) Contopoulos, G., Kandrup, H.E., Kaufmann, D.: Fractal properties of escape from a twodimensional potential. Phys. D 64, 310323 (1993)
 (22) Contopoulos, G., Harsoula, M., LukesGerakopoulos, G.: Periodic orbits and escapes in dynamical systems. Celest. Mech. Dyn. Astron. 113, 255278 (2012)
 (23) de Moura, A.P.S., Letelier, P.S.: Fractal basins in HÃ©nonHeiles and other polynomial potentials. Phys. Lett. A 256, 362368 (1999)
 (24) Deprit, A.: The Lissajous transformation. I. Basics. Celest. Mech. Dyn. Astron. 51, 202225 (1991)
 (25) Deprit, A., Elipe, A.: The Lissajous transformation. II. Normalization. Celest. Mech. Dyn. Astron. 51, 227250 (1991)
 (26) Eckhardt, B.: Irregular scattering. Phys. D 33, 8998 (1988)
 (27) Elipe, A.: Complete reduction of oscillators in resonance . Phys. Rev. E 61, 64776484 (2000)
 (28) Elipe, A., Deprit, A.: Oscillators in resonance. Mech. Res. Commun. 26, 635640 (1999)
 (29) Ernst, A., Just, A., Spurzem, R., Porth, O.: Escape from the vicinity of fractal basin boundaries of a star cluster. Mon. Not. R. Astron. Soc. 383, 897906 (2008)
 (30) Ferrer, S., Lara, M., Palacián, J., Juan, J.S., Viartola, A., Yanguas, P.: The HénonHeiles problem in three dimensions. I. Periodic orbits near the origin. Int. J. Bifurc. Chaos 8, 11991213 (1998)
 (31) Ferrer, S., Lara, M., Palacián, J., Juan, J.S., Viartola, A., Yanguas, P.: The HénonHeiles problem in three dimensions. II. Relative equilibria and bifurcations in the reduced problem. Int. J. Bifurc. Chaos 8, 12151229 (1998)
 (32) Giorgilli, A., Galgani, L.: From integrals from an autonomous Hamiltonian system near an equilibrium point. Celest. Mech. 17, 267280 (1978)
 (33) Hénon, M., Heiles, C.: The applicability of the third integral of motion: Some numerical experiments. Astron. J. 69, 7379 (1964)
 (34) Jung, C., Scholz, H.: Cantor set structures in the singularities of classical potential scattering. J. Phys. A 21, 36073617 (1988)
 (35) Kandrup, H.E., Siopis, C., Contopoulos, G., Dvorak, R.: Diffusion and scaling in escapes from twodegreesoffreedom Hamiltonian systems. Chaos 9, 381392 (1999)
 (36) Karanis, G.I., Vozikis, Ch.L.: Fast detection of chaotic behavior in galactic potentials. Astron. Nachr. 329, 403412 (2007)
 (37) Kennedy, J., Yorke, J.A.: Basins of Wada. Phys. D 51, 213225 (1991)
 (38) Motter, A.E., Lai, Y.C.: Dissipative chaotic scattering. Phys. Rev. E 65 R0152051–4 (2002)
 (39) Ott, E., Tél, T.: Chaotic scattering: An introduction. Chaos 3, 417426 (1993)
 (40) Petit, J.M., Hénon, M.: Satellite encounters. Icarus 66, 536555 (1986)
 (41) Poon, L., Campos, J., Ott, E., Grebogi, C.: Wada basins boundaries in chaotic scattering. Int. J. Bifurc. Chaos 6, 251266 (1996)
 (42) Press H.P., Teukolsky S.A, Vetterling W.T., Flannery B.P.: Numerical Recipes in FORTRAN 77, 2nd Ed., Cambridge Univ. Press, Cambridge, USA (1992)
 (43) Saito, N., Ichimura, A.: In: Casati, G., Ford, J. (eds.) Stochastic Behavior in Classical and Quantum Hamiltonian Systems. Lecture Notes in Physics, vol. 93, p. 137. Springer, Berlin (1979)
 (44) Seoane, J.M.: Effects of weak perturbations in open Hamiltonian systems, PhD thesis. Universidad Rey Juan Carlos, Spain (2007)
 (45) Seoane, J.M., Aguirre, J., Sanjuán, M.A.F., Lai, Y.C.: Basin topology in dissipative chaotic scattering. Chaos 16, 0231011–8 (2006)
 (46) Seoane, J.M., Sanjuán, M.A.F., Lai, Y.C.: Fractal dimension in dissipative chaotic scattering. Phys. Rev. E 76, 0162081–6 (2007)
 (47) Seoane, J.M., Sanjuán, M.A.F.: Exponential decay and scaling laws in noisy chaotic scattering. Phys. Lett. A 372, 110116 (2008)
 (48) Seoane, J.M., Huang, L., Sanjuán, M.A.F., Lai, Y.C.: Effects of noise on chaotic scattering. Phys. Rev. E 79, 0472021–4 (2009)
 (49) Seoane, J.M., Sanjuán, M.A.F.: Escaping dynamics in the presence of dissipation and noisy in scattering systems. Int. J. Bifurc. Chaos 9, 27832793 (2010)
 (50) Siopis, C.V., Contopoulos, G., Kandrup, H.E.: Escape probabilities in a Hamiltonian with two channels of escape. New York Acad. Sci. Ann. 751, 205212 (1995)
 (51) Siopis, C.V., Kandrup, H.E., Contopoulos, G., Dvorak, R.: Universal properties of escape. New York Acad. Sci. Ann. 773, 221230 (1995)
 (52) Siopis, C.V., Kandrup, H.E., Contopoulos, G., Dvorak, R.: Universal properties of escape in dynamical systems. Celest. Mech. Dyn. Astron. 65, 57681 (1996)
 (53) Skokos, C.: Alignment indices: A new, simple method for determining the ordered or chaotic nature of orbits. J. Phys. A: Math. Gen. 34, 1002910043 (2001)
 (54) Skokos, Ch., Antonopoulos, Ch., Bountis, T.C., Vrahatis, M.N.: Detecting order and chaos in Hamiltonian systems by the SALI method. J. Phys. A: Math. Gen. 37, 62696284 (2004)
 (55) Zotos, E.E.: Application of new dynamical spectra of orbits in Hamiltonian systems. Nonlin. Dyn. 69, 20412063, (2012)
 (56) Zotos, E.E.: Trapped and escaping orbits in axially symmetric galactictype potential. PASA 29, 161â173 (2012)
 (57) Zotos, E.E.: The Fast Norm Vector Indicator (FNVI) method: A new dynamical parameter for detecting order and chaos in Hamiltonian systems. Nonlin. Dyn. 70, 951978 (2012)
 (58) Zotos, E.E.: Revealing the evolution, the stability and the escapes of families of resonant periodic orbits in Hamiltonian systems. Nonli. Dyn. 73, 931962 (2013)
 (59) Zotos, E.E., Caranicolas, N.D.: Are seminumerical methods an effective tool for locating periodic orbits in 3D potentials? Nonlin. Dyn. 70, 279287 (2012)
 (60) Zotos, E.E., Caranicolas, N.D.: Order and chaos in a new 3D dynamical model describing motion in nonaxially symmetric galaxies. Nonlin. Dyn. 74, 12031221 (2013)