Methods in chaos detection

# Chaos detection tools: application to a self–consistent triaxial model

## Abstract

Together with the variational indicators of chaos, the spectral analysis methods have also achieved great popularity in the field of chaos detection. The former are based on the concept of local exponential divergence. The latter are based on the numerical analysis of some particular quantities of a single orbit, e.g. its frequency. In spite of having totally different conceptual bases, they are used for the very same goals such as, for instance, separating the chaotic and the regular component. In fact, we show herein that the variational indicators serve to distinguish both components of a Hamiltonian system in a more reliable fashion than a spectral analysis method does.

We study two start spaces for different energy levels of a self–consistent triaxial stellar dynamical model by means of some selected variational indicators and a spectral analysis method. In order to select the appropriate tools for this paper, we extend previous studies where we make a comparison of several variational indicators on different scenarios. Herein, we compare the Average Power Law Exponent (APLE) and an alternative quantity given by the Mean Exponential Growth factor of Neary Orbits (MEGNO): the MEGNO’s Slope Estimation of the largest Lyapunov Characteristic Exponent (SElLCE). The spectral analysis method selected for the investigation is the Frequency Modified Fourier Transform (FMFT).

Besides a comparative study of the APLE, the Fast Lyapunov Indicator (FLI), the Orthogonal Fast Lyapunov Indicator (OFLI) and the MEGNO/SElLCE, we show that the SElLCE could be an appropriate alternative to the MEGNO when studying large samples of initial conditions. The SElLCE separates the chaotic and the regular components reliably and identifies the different levels of chaoticity. We show that the FMFT is not as reliable as the SElLCE to describe clearly the chaotic domains in the experiments. We use the latter indicator as the main variational indicator to analyse the phase space portraits of the model under study.

###### keywords:
methods: numerical – galaxies: kinematics and dynamics.
45

## 1 Introduction

The understanding of a realistic dynamical model is strongly related to the capability of identifying the chaotic and the regular nature of its orbits.

One the one hand, one of the main features of chaotic orbits is their high sensitivity to the initial conditions (hereinafter i.c.), and the concept of local exponential divergence applies perfectly to identify such correspondence. On the other hand, the regular orbits do not show this behaviour. A seminal contribution to the field of chaos detection has been made by Lyapunov when he introduced the Lyapunov Characteristic Exponents, hereinafter LCEs (Oseledec, 1968; Lyapunov, 1992). The LCEs are theoretical quantities that measure the local rate of exponential divergence. Therefore, the LCEs are appropriate to distinguish between chaotic and regular motion.

The LCEs constitute the backbone of variational indicators of chaos (VICs). Moreover, their numerical implementation (Benettin et al., 1980; Froeschlé, 1984; Tancredi et al., 2001; Skokos, 2010), i.e. the Lyapunov Indicators (LIs), was a major contribution to the development of both fast and easy–to–compute new VICs. Nowadays, there is a plethora of VICs in the literature, such as the Smaller Alignment Index, SALI (Skokos, 2001; Skokos et al., 2004; Széll et al., 2004; Bountis & Skokos, 2006; Carpintero, 2008; Antonopoulos et al., 2010), and its generalized version, the Generalized Alignment Index, GALI (Skokos et al., 2007, 2008; Manos & Athanassoula, 2011); the Mean Exponential Growth factor of Nearby Orbits, MEGNO (Cincotta & Simó, 2000; Cincotta et al., 2003; Giordano & Cincotta, 2004; Goździewski et al., 2005; Gayon & Bois, 2008; Lemaître et al., 2009; Hince et al., 2010; Compère et al., 2011; Maffione et al., 2011a); the Fast Lyapunov Indicator, FLI (Froeschlé et al., 1997a, b; Froeschlé & Lega, 1998, 2000; Lega & Froeschlé, 2001; Guzzo et al., 2002; Froeschlé & Lega, 2006; Paleari et al., 2008; Todorović et al., 2008; Lega et al., 2010); its variant, the Orthogonal Fast Lyapunov Indicator, OFLI (Fouchard et al., 2002); and a second order variant of the OFLI, the OFLI (Barrio, 2005; Barrio et al., 2009a, b, 2010). Furthermore, we can include the VICs based on the properties of dynamical spectra like the invariant spectra of Stretching Numbers or Local Lyapunov Characteristic Numbers, LLCNs (Contopoulos et al., 1997; Lega & Froeschlé, 1998), the invariant spectra of Helicity or Twist Angles (Voglis & Contopoulos, 1994; Contopoulos & Voglis, 1996, 1997; Contopoulos et al., 1997; Froeschlé & Lega, 1998; Lega & Froeschlé, 1998; Voglis et al., 1998) and the Spectral Distance (Voglis et al., 1999). Finally, we include the Relative Lyapunov Indicator, RLI (Sándor et al., 2000, 2004, 2007; Széll et al., 2004), which is not based on the evolution of the solution of the first variational equations as the others, but on the evolution of two different but very close orbits.

Notwithstanding the large number of VICs cited, we could also include alternative methods, for instance, the Average Power Law Exponent, APLE (Lukes-Gerakopoulos et al., 2008), which is based on the concept of Tsallis entropy; evolutionary algorithms (Petalas et al., 2009); and the record could indeed continue.

Together with the VICs, the spectral analysis methods (SAMs), based on the numerical analysis of some particular quantities of a single orbit, e.g. its frequency, have also achieved great popularity in the field of chaos detection. Among the SAMs we find the method outlined by Laskar (1990): the Frequency Map Analysis, FMA (Laskar et al., 1992; Laskar, 1993; Papaphilippou & Laskar, 1996, 1998) and a modification of the latter in order to improve the precision in the computed frequencies and amplitudes, the Frequency Modified Fourier Transform, FMFT (Sidlichovský & Nesvorný, 1997) or other alternatives using Fast Fourier Transform (FFT) techniques (Michtchenko & Ferraz–Mello, 1995; Ferraz–Mello et al., 2005; Michtchenko et al., 2010), among others. As an easy identification of the important resonances (regions of phase space occupied by orbits with commensurable orbital frequencies) is possible with the SAMs, they are able to identify the major resonant orbit families. The relative importance of each family regarding the phase space can be assessed from the number of orbits associated with the particular resonance the family belongs to. That is, the SAMs are not only suitable to determine the resonance web like the VICs (Kaneko & Konishi, 1994; Cincotta et al., 2003; Froeschlé et al., 2006; Lukes-Gerakopoulos et al., 2008) but also to identify the individual resonances by their frequency vectors.

Though the conceptual bases of the VICs and the SAMs are totally different, many researchers use both types of techniques for similar goals. Nevertheless, we show herein that this seems to be a no reliable approach to study some important aspects of the dynamics of a Hamiltonian system. Moreover, in the literature there are plenty of evidence that the complementary use of both types of techniques results an efficient and solid way to gather dynamical information from a Hamiltonian system (Froeschlé et al., 1997b; Froeschlé & Lega, 1998; Guzzo et al., 2002; Kalapotharakos & Voglis, 2005; Barrio et al., 2009b). For instance, in Froeschlé et al. (1997b) the authors compare the sensitivity of the FLI and the FMA on the standard map and on the Hénon & Heiles (1964) potential. In Guzzo et al. (2002) the authors use the FLI and the Analytically Filtered Fourier Analysis, AFFA (Guzzo & Benettin, 2001), another example of SAM, to determine the value of the critical parameter at which the transition from the unfortunately called Nekhoroshev to the Chirikov regime occurs in a quasi–integrable Hamiltonian model and a standard –dimensional map. In Kalapotharakos & Voglis (2005), the authors use a measure related with the LCEs, the alignment indices (Voglis et al., 1998, 1999; Skokos, 2001) and an index computed with the FMFT to study two different N–body models simulating elliptical galaxies.

One of the aspects that makes a research accurate is the reliability of the tools available for the analysis. This reliability is a measure of, for instance, the confidence and efficiency of such tools. Regarding chaos detection, we believe that it is on the diversity of the employed techniques that the fruitfulness of the research relies. Therefore, we started a series composed of four papers in order to achieve the best combination of tools to study a given dynamical problem. In a first paper, Maffione et al. (2011a), we test the MEGNO against the most widely used chaos detection technique in the literature, which is the LI, on a somewhat realistic problem: a self–consistent triaxial stellar (ScTS) dynamical model (Muzzio et al., 2005; Cincotta et al., 2008). In the second and third papers, Maffione et al. (2011b) (hereinafter M11) and Darriba et al. (2012) (hereinafter D12), we make a comparison of most of the VICs mentioned earlier on rather simple problems, i.e. on mappings and on the model of Hénon & Heiles (1964), respectively. Our purpose here is to extend the previous comparisons of VICs to alternative techniques, including a spectral analysis method. Further, this last paper of the series is regarded as a closing report on the subject and the conclusions comprise the whole investigation. We continue with the ScTS model and study two start spaces for four different energy levels by means of several chaos detection tools. In order to select appropriate VICs for the study, we consider the results of M11 and D12 and add a rather new technique, the APLE, and an alternative quantity given by the MEGNO: the MEGNO’s Slope Estimation of the largest LCE (SElLCE) to their comparison. The MEGNO is part of the suggested VICs to compose the CIsF6 in D12. In fact, we are interested in studying if its variant, i.e. the SElLCE, is more reliable than the MEGNO on its own.

The study is organized as follows: in Section 2 we describe the ScTS model and introduce the and the start spaces. In Section 3, we briefly introduce the chaos detection techniques to be used within the reported research. We start with the VICs: the LI, the MEGNO and the SElLCE (Section 3.1); the FLI and the OFLI (Section 3.2) and the APLE (Section 3.3). We finish with a short FMFT introduction (Section 3.4). In Section 4, we deal with the capability of the FMFT to resolve the different structures in a divided phase space (Section 4.1) and compare its performance with the SElLCE in both start spaces of the ScTS model (Section 4.2). In Section 5, we study the efficiency of several VICs on the start space to decide which ones are appropriate to investigate the dynamics of the ScTS model (Sections 5.1 and 5.2). In Section 6 we analyse for four different energy surfaces, the (Section 6.1) and the (Section 6.2) start spaces by means of a complementary use of the VICs selected in the previous section. Finally, Section 7 summarizes the main conclusions of the whole investigation.

As the extensive used of acronyms in the text might be overwhelming for the reader, in Table 1 we list in alphabetical order the most important acronyms used in this paper.

## 2 The model and the start spaces

In order to obtain a fairly realistic dynamical model of an elliptical galaxy, in Muzzio et al. (2005) the authors follow the cold disipationaless collapse of 100000 particles randomly distributed following a law within a sphere. The velocities were randomly chosen from a spherical Gaussian distribution but only their tangential component was retained. After the system had relaxed, there remained 86818 particles resembling an elliptical galaxy: the ScTS model. The system observes the de Vaucouleurs law shown by Muzzio et al. (2005) in their Fig. 2. The model reproduces many dynamical characteristics of real elliptical galaxies, such as mass distribution, flattening, triaxiality and rotation (Muzzio, 2006). The ScTS model has a strong triaxiality and a flattening that increases from the border of the system to its center (see Table I in the same paper). The resulting triaxial potential has semi–axis satisfying the condition , and its minimum, which is close to , matches the origin. The potential is less flattened than the mass distribution, as expected. See Table I in Muzzio et al. (2005).

The ScTS model is an –body potential which is neither smooth nor stationary. As the authors needed to compute the orbits and their corresponding LIs, they froze the potential and represented it with a quadrupolar approximation. The equation that reproduces the potential is:

 V(x)=−f0(x)−fx(x)⋅(x2−y2)−fz(x)⋅(z2−y2)

with and where

 fn(x)=αn[pann+δann]acnan, (1)

with the square of the softened radius given by when , or for , and , , , are constants.

The selected value for the softening parameter is for any . The adopted values for the constants , , and are given in Table 2.

The stationary character of the parameters given in Table 2 were tested by performing several fits at different times after virialization, resulting in a precision of . For further details on the ScTS model refer to Muzzio et al. (2005); Cincotta et al. (2008).

We have already used the same model to test the MEGNO in a previous study (Maffione et al., 2011a).

The spaces of i.c. used in this investigation are the and the start spaces for four energy surfaces: and . The start space is a set of i.c. of the form with while is restrained by the energy condition (Papaphilippou & Laskar, 1998). The start space is a set of i.c. of the form with and is restrained by the energy condition (Schwarzschild, 1993). On the one hand, the combination of both start spaces is an adequate choice to sample many of the different orbital families within the triaxial model (Schwarzschild, 1993; Papaphilippou & Laskar, 1998). On the other hand, the stable boxlets (i.e. resonant box orbits) which are centrophobic are not well represented in our start space.

Therefore, the and the start spaces of the ScTS model seem to provide fairly realistic scenarios for testing many characteristics of both the VICs of the package and the FMFT.

## 3 The techniques

The aim of this section is to briefly introduce the techniques we will use in the forthcoming study.

### 3.1 The LI, the MEGNO and the SElLce

Considering a continuous dynamical system defined on a differentiable manifold , where characterises the state of the system at time , being the state of the system at time , . The state of the system after two consecutive time steps and will be given by the composition law: .

Moreover, the tangent space of maps onto the tangent space of according to the operator and following the rule , where is an initial deviation vector (hereinafter i.d.v.). The action of such operator at consecutive time intervals satisfies the equation:

 dxΦt+t′=dΦt′(x)Φt∘dxΦt′.

If we suppose that our manifold has some norm denoted by , we can define the useful quantity:

 λt(x)=∥dxΦt(w)∥∥w∥

called “growth factor” in the direction of .

Let with be an –dimensional Hamiltonian, that we suppose autonomous just for the sake of simplicity. Introducing the following notation:

 x=(p,q)∈R2N,f(x)=(−∂H/∂q, ∂H/∂p)∈R2N,

the equations of motion can be written in a simple way as

 ˙x=f(x). (2)

Let be an arc of an orbit of the flow (2) over a compact energy surface: , , then

 γ(x0;t)={x(x0;t′):x0∈Mh, 0≤t′

We can gain fundamental information about the Hamiltonian flow in the neighborhood of any orbit through the largest LCE: lLCE, defined as:

 χγ=χ[γ(x0;t)]=limt→∞1tlnλt[γ(x0;t)] (3)

with

 λt[γ(x0;t)]=∥dγΦt(w)∥∥w∥

where is an “infinitesimal displacement” from at time . The fact that the lLCE (and its truncated value, the so–called LI for T finite) measures the local mean exponential rate of divergence of nearby orbits is clearly understood when Eq. (3) is written in an integral fashion:

 χγ=limt→∞1t∫t0∥˙dγΦt′(w)∥∥dγΦt′(w)∥dt′=¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯(∥˙dγΦt′(w)∥∥dγΦt(w)∥), (4)

where the bar denotes time average.

Finally, the orbit is classified as regular or chaotic if the LI tends to zero (as ) or to a positive value, respectively (a detailed discussion on the theory and computation of the LCEs can be found in the review of Skokos (2010)).

Now, we can introduce the quantity through the expression:

 Yγ=Y[γ(x0;t)]=2t∫t0∥˙dγΦt′(w)∥∥dγΦt′(w)∥t′dt′,

which is related to the integral in Eq. (4); i.e., in case of an exponential increase of , , the quantity can be considered as a weighted variant of the integral in Eq. (4). Notice that the quantity verifies for regular motion and for chaotic motion, with , which show that, in case of regular motion converges to 0 faster than does (which goes to zero as , while for chaotic motion both magnitudes approach the positive lLCE at a rather similar rate.

Instead of using the instantaneous rate of increase, , we average the logarithm of the growth factor, . Introducing the time average (MEGNO):

 ¯¯¯¯Yγ=¯¯¯¯Y[γ(x0;t)]≡1t∫t0Yγdt′,

we notice that the time evolution of the MEGNO can be briefly described in a suitable and unique expression for all kinds of motion. Indeed, its asymptotic behaviour can be summarized as follows:

 ¯¯¯¯Yγ≈aγt+bγ, (5)

where and for quasi–periodic motion, while and for irregular, stochastic motion. Deviations from the value (and ) indicate that is close to some particular objects in phase space, being or for stable periodic orbits (resonant elliptic tori), or unstable periodic orbits (hyperbolic tori), respectively.

Refer to Fig. 1 (c) and (d) from Cincotta et al. (2003) to see the general behaviour of the MEGNO for regular and chaotic orbits.

It is possible to estimate the LI of the orbit from Eq. (5) if we apply a simple linear least–squares fitting on the MEGNO (Cincotta & Simó, 2000; Cincotta et al., 2003). This estimation is the earlier mentioned SElLCE indicator. The least–squares fitting at the end of the integration process uses the last 80% of the points in order to avoid the initial transient. Thus, the SElLCE almost considers the full history of the orbit which makes it a very sensitive VIC.

The performances of the MEGNO’s time evolution curves to characterise the different levels of stability and chaoticity of the orbits shown in M11 and D12 was excellent. On the other hand, the description of the regular components in a divided phase space by means of the final values of the MEGNO (i.e. the values of the indicator at the end of the total integration time), was not satisfactory. The SElLCE seems to be an appropriate alternative to solve the problem.

Fig. 1 shows the performances of the LI, the MEGNO and the SElLCE for two sets of i.c., one of them located in the start space (top panel) and the other set located in the start space (bottom panel). In order to ensure that the VICs are well computed for a given orbit, the total integration time should verify , where the function is some characteristic time–scale which depends on the energy surface . The function is taken as the period of the stable X–axis periodic orbit and we consider the energy surface for both sets of i.c. ( units of time –hereinafter u.t.– see Maffione et al. (2011a)). Therefrom, on fixing the condition to obtain reliable values for the VICs, we conclude that a total integration time of u.t. would be appropriate for the experiment (Maffione et al., 2011a). In the start space, the set is composed of 790 i.c. of the form and . In the start space, the set is composed of 646 i.c. of the form and .

In Fig. 1, the SElLCE and the LI show very similar final values for chaotic orbits. The MEGNO grows exponentially fast and reaches higher final values () than the SElLCE or the LI. The MEGNO clearly distinguishes between different levels of chaoticity. Consequently, the SElLCE seems to be an unnecessary alternative to the MEGNO for describing the chaotic component. In case of a divided phase space, where regular and chaotic components are mixed, the asymptotically behaviour of MEGNO’s final values for regular orbits may hide their different stability levels (see M11 and D12). The SElLCE tends to zero for regular orbits. Moreover, the SElLCE tends to decrease faster than the LI for regular motion (Cincotta & Simó, 2000; Cincotta et al., 2003). So, no matter if the phase space is strongly divided, this behaviour allows the SElLCE’s final values to show certainly (and faster than the LI) different levels of stability for the regular orbits. Thus, it seems to be appropriate to use the time evolution curves of the MEGNO for the individual analysis of the orbits (M11 and D12) and the final values of the SElLCE as a global VIC to obtain the portraits of divided phase spaces (see Sections 4.2, 6.1 and 6.2).

The MEGNO has a low computational cost and the linear least–squares fitting in the computation of the SElLCE does not require a significant amount of extra time. Hence, the SElLCE is an alternative quantity that has also a low computational cost.

### 3.2 The FLI and the OFLI

The FLI is a quantity intimately related to the lLCE (Froeschlé et al., 1997b, a; Froeschlé & Lega, 2000; Lega & Froeschlé, 2001; Guzzo et al., 2002; Froeschlé & Lega, 2006) and the MEGNO (Mestre et al., 2011). The FLI is able to distinguish between regular and chaotic (weakly chaotic) motion (Froeschlé et al., 1997b, a) and also between resonant and non–resonant motion (Froeschlé & Lega, 2000; Lega & Froeschlé, 2001; Guzzo et al., 2002) using only the first part of the computation of the lLCE.

Here, we use the initial definition of the FLI given in Froeschlé et al. (1997b, a) like in D12, i.e. we use a basis of i.d.v. for its computation. The definition given in Froeschlé & Lega (2000) uses only one i.d.v., the results do not change. Nevertheless, the computing time is obviously reduced (Froeschlé & Lega, 2000).

On an –dimensional Hamiltonian , we follow the time evolution of unit i.d.v. (e.g. the ScTS model is –dimensional and the basis consists of 6 i.d.v.). The FLI at time is defined by the highest norm among the unit i.d.v. of the basis, as follows:

 FLI(t)=supt[∥w(t)1∥,∥w(t)2∥,…,∥w(t)2N∥].

For further details, refer to Froeschlé et al. (1997b, a).

For both regular and chaotic motion, the FLI tends to infinity as time increases, but with completely different rates. The FLI grows linearly with time for regular motion and it grows exponentially fast for chaotic motion. Furthermore, the FLI grows linearly with time, but with different rates in the case of resonant regular orbits and non–resonant ones. Although the curves corresponding to resonant and non resonant motions are separated, the oscillations due to the distortion of the orbits themselves may hide the distinction between the libration islands and the tori. Therefore, in Froeschlé & Lega (2000) the authors not only introduce the FLI using only one i.d.v. but also replace the definition of the FLI by its running average. We also applied an average to smooth the curves of the FLI.

In the case of the OFLI7, we take the orthogonal component for the flow for each unit i.d.v. () of the basis at every time step and retain the largest component among them:

 OFLI(t)=supt[w(t)⊥1,w(t)⊥2,…,w(t)⊥2N].

This modification makes the OFLI a VIC which can easily distinguish periodicity among the orbits of the regular component. The OFLI for periodic orbits oscillates around a constant value, while for quasiperiodic motion and chaotic motion it has the same behaviour as the FLI.

Refer to Figs. 5 and 7 from Fouchard et al. (2002) to see the general behaviour of the FLI and the OFLI for regular and chaotic orbits.

### 3.3 The APLE

In the early study of Tsallis et al. (1997) the authors show that when a nonlinear dynamical system is in the regime of the so–called edge of chaos, the q–entropy rate of increase remains constant for a long time (metastable states). Moreover, they discuss that this behaviour may be associated with a power–law rather than with exponential sensitivity of the orbits to the i.c. The interested reader may refer to Tsallis et al. (2002) for a “pedagogical guided tour” on the subject.

In Lukes-Gerakopoulos et al. (2008) the authors find that the time evolution of the i.d.v. for a weakly chaotic orbit (i.e. a combination of a linear and an exponential law: , with and ) justifies theoretically why metastable states with a constant rate of increase of the q–entropy appear in the nonlinear dynamical system. They argue that this law is reflected upon a transient behaviour in which the growth of the i.d.v. is almost a power law: , with (and the exponent can be associated to a –exponent via the relationship ). In fact, they introduce a method to compute the –exponent and use it as a VIC (the APLE) in order to distinguish regular orbits from nearby weakly chaotic orbits.

For an –dimensional Hamiltonian consider a partitioning of the –dimensional phase space into a large number of volume elements of size for some small and let be the initial condition of an orbit located in a particular volume element. The authors introduce the APLE: , as follows:

 pγ=ln(|w(t)|2|w(t1)|2)2ln(tt1),

where and is one of the deviation vectors of an orthogonal basis of the tangent space to at the initial point . Every has length greater or equal to and is a transient initial time of orbit evolution.

Therefore, in case of regular motion, tends to 1, while for chaotic motion it grows exponentially fast. Moreover, the APLE shows oscillations around a constant value for a transient time interval before the exponential growth for weakly chaotic orbits. The length of the time interval is directly connected with the hyperbolicity of the local phase space.

Refer to Fig. 2 from Lukes-Gerakopoulos et al. (2008) to see the general behaviour of the APLE for regular and chaotic orbits.

### 3.4 The FMFT

The detailed description of the FMFT8 has been published on Sidlichovský & Nesvorný (1997). However, we should mention the conceptual aspects of this technique. Since the main advantage of the FMFT is that it improves the computation of the amplitudes and the frequencies given by the FMA (Laskar, 1990; Laskar et al., 1992; Laskar, 1993, 1996; Papaphilippou & Laskar, 1996, 1998; Laskar, 2003), we decided to introduce only a short description of the latter.

Let be a non–integrable Hamiltonian which is also non–degenerate:

 det(∂\mn@boldsymbolν(I)∂I)=det(∂2H(I)∂I2)≠0,

with the action-like variables and the frequency vector. The KAM theory states that for sufficient small values of the perturbation parameter , there exists a Cantor set of that satisfies the Diophantine condition:

 ∥k⋅\mn@boldsymbolν∥>κϵ|k|m,

where and are constants and is the norm of . Furthermore, the system still possesses invariant tori with linear flow (KAM tori). The solutions lie on these invariant tori and are given in complex variables through their Fourier series:

 zj(t)=zj0eiνj×t+∑mam(\mn@boldsymbolν)ei⟨m,\mn@boldsymbolν⟩

where the coefficients depend smoothly on the frequencies . If we keep all the angle–like variables fixed (, is the –dimensional torus, fixed to ), we obtain a frequency map in defined as:

 F\mn@boldsymbolθ0:B→Ωϵ;I→p2(Ψ−1(\mn@boldsymbolθ0,I)) (6)

where is the projection on . The FMA numerically gives through an iterative process the frequency map (i.e. amplitudes and frequencies) defined over the whole domain , which coincides, within the numerical precision adopted, with the of Eq. (6) in the set of KAM tori. The frequency map is obtained seeking for the quasiperiodic approximations of the solutions over a finite time–span through a finite series of terms:

 zj(t)=zj0eiνj×t+s∑k=1amkei⟨mk,\mn@boldsymbolν⟩. (7)

Once we obtain the quasiperiodic approximation from Eq. (7), we are able to establish a correspondence between the action–like variables () and the rotation numbers (or frequency vector ). This correspondence is the so–called FMA (Laskar, 1990). Regular regions of the phase space allow a very accurate estimation of real rotation numbers in a given time interval (for instance, using the MFT with the Hanning window, using the MFT without the Hanning window, using the usual FFT techniques, where is set as the whole time–span), and should not differ from the rotation numbers obtained in another time interval, considering a certain accuracy. Strictly speaking, the frequencies are only defined on KAM tori, the FMA numerically estimates over a finite time–span a frequency vector for any initial condition, though. If the estimates of the rotation numbers given by the FMA vary (above the required accuracy) between different time periods, it implies that the corresponding KAM tori are destroyed. Thus, the study of this frequency map (i.e. the FMA) regarding its constancy in time provides important clues about the dynamics of the system.

For a comprenhensive discussion of the FMA refer to Laskar (2003).

In order to calculate the frequencies with the FMFT, we need to compute the equations of motion in the ScTS model. We use the Taylor method (Jorba & Zou, 2005) that proved to be very convenient for the task in D12. The precision in the computation of the coordinates is . The numerical computation of the VICs was done using the DOPRI8 routine (Prince & Dormand, 1981). This routine achieved better results than Taylor when we required the simultaneous integration of the coupled system of equations of motion and variational equations in the ScTS model (see also D12). The conservation of the energy with the DOPRI8 routine was .

All the computations in the forthcoming investigation were done using the following configuration. Hardware: CPU, 2 x Dual XEON 5450, Dual Core 3.00GHz; M.B., Intel S5000VSA; RAM, 4GB(4x1GB), Kingston DDR–2, 667MHz, Dual Channel. Software: gfortran 4.2.3.

## 4 The FMFT as a global indicator of chaos

In this section, we compare the performances of the FMFT and the SElLCE as global indicators of chaos.

### 4.1 Accuracy in the estimation of the frequencies

In order to analyse the nature of each orbit we apply the FMFT to the following function in each degree of freedom (hereinafter, d.o.f.):

 ψj(tk)=ψRj(tk)+iψCj(tk) (8)

with () a complex function, where (the real part of ) corresponds to the Cartesian coordinates ( and ) and (the conjugate part of ) corresponds to the Cartesian velocities ( and ). The with are the trajectory samplings. Every set of points generates a quasiperiodic approximation and we keep the frequency with the maximum amplitude. We call this frequency the fundamental frequency of the system in the d.o.f.

Theoretically speaking, if the fundamental frequencies correspond to a regular orbit, they should remain constant in time. Numerically speaking, the variation exists and it is due to the precision of the computation. Thus, in order to determine such a precision for our particular experiment, we take samples of regular orbits, according to both convergent values of the LI and the MEGNO ( or , respectively, see Maffione et al. (2011a)).

We are aware that samples obtained from the use of VICs can only be statistically rich in regular orbits. That is, the actual values of the VICs are reached only at infinite time and as we have just truncated approximations of those values, some of the orbits might be missclassified. In order to be as certain as possible about the dominant type of orbits from the samples, we follow their behaviour for (see Section 3.1). In a galaxy like the one represented by the ScTS model, the crossing time is u.t. If we suppose that the Hubble time is of the order of crossing times (a high value indeed), then we are integrating the orbits for and Hubble times for the energy surfaces ( u.t.) and , respectively. Therefore, seems to be a lapse of time large enough to provide a physical meaningful characterization of the motion of the orbits and reliable convergent values of the VICs for most of the orbits in the sample. Finally, for all practical purposes, most of the orbits of the samples are regular orbits.

The samples are taken in the and the start spaces for four different energy surfaces, namely: , ( u.t.), ( u.t.) and (i.e. a total of eight samples). There are 100 i.c. in each sample and we integrate the equations of motion for for every i.c. Although is the time–span we use to compute the VICs reliably, the results in the determination of the fundamental frequencies with the FMFT are not very different using only . Besides, the CPU time is certainly reduced. Therefore, we apply the FMFT on two overlapping time intervals of (Wachlin & Ferraz–Mello, 1998)9. The number of points used in the estimation of the frequencies is which makes points for each . Finally, we estimate the fundamental frequencies for each time interval and calculate the differences (the accuracy in the determination of the frequencies should be according to Section 3.4 and Muzzio (2006)). The concomitant results for the energy surfaces and for both the and the start spaces are shown in Fig. 2.

In Fig. 2 we show that the differences in the computed values of the fundamental frequencies between both overlapping time intervals for regular orbits (actually, some of them might be chaotic) is not generally below the theoretical value, particularly in the start space. Thus, there exists a range of values where the regular orbits could not be reliably distinguished from the chaotic orbits by means of the computation of the fundamental frequencies with the FMFT.

In Table 3 we show the corresponding arithmetic means and standard deviations of the variations of the fundamental frequencies for the energy surfaces: , , and on the start space. The arithmetic means show an average conservation between the 7th and the 8th decimal in the fundamental frequencies. Nevertheless, the standard deviations oscillate around an average of the 5th decimal. Therefore, there is an interval where the FMFT can not distinguish clearly between regular and chaotic orbits in the experiment. However, we can infer certain limits which can help us to identify if the orbit is regular or chaotic using the FMFT. For instance, the variations in the computed fundamental frequencies for the regular orbits only affect decimals higher than the 4th (see Fig. 2). Thus, if there is a difference in the computation of the fundamental frequencies of the orbit in decimals lower than the 4th, it is highly probable that the orbit is chaotic. This result is in complete agreement with the one shown in Muzzio (2006) for the same ScTS model and with the practical limit used also in Aquilano et al. (2007); Muzzio et al. (2009).

The results for the start space are similar than those for the start space. Nevertheless, for lower energy values (Fig. 2, bottom right panel), the averaged precision in the computation of the frequencies improved substantially.

We also make some experiments with the FMFT, using synthetic signals for verifying the precision in the computation of the fundamental frequencies. The results were remarkable better and are closer to the theoretical estimates (Muzzio, 2006). Nevertheless, this is not the case for the relatively complex ScTS model and the precision in the estimations is clearly diminished. Thus, we use the same model and the same spaces we study to calibrate the FMFT.

### 4.2 A comparative evaluation between the FMFT and the SElLce

We are interested in testing the performance of the FMFT as a global indicator of chaos. Thus, we apply it together with a variational indicator in order to characterise the phase space portraits of two different regions. One of the regions is located in the start space (the region is identified in the top panel of Fig. 3), whereas the other is located in the start space (which is identified in the bottom panel of the same figure). Both regions belong to the energy surface . The VIC selected for this study is the SElLCE (see Section 3.1).

In order to apply the FMFT as a global indicator of chaos, we follow the same procedure as described in Section 4.1. Thus, we consider two overlapping time intervals of periods. The first time interval is u.t. and the second time interval is u.t. According to Section 3.4, the accuracy in the determination of the frequencies should be . The computation is over 624100 i.c., within the start space, covering the region depicted on the top panel of Fig. 3 and over 596258 i.c., within the start space, covering the region depicted on the bottom panel of the same figure.

To use the FMFT as a global indicator of chaos and compute the corresponding phase space portraits, we need to define a new quantity (Wachlin & Ferraz–Mello, 1998). Such a quantity is , being , where is the frequency (computed with the FMFT) associated with the d.o.f “j” for the time interval “(i)”. All the fundamental frequencies have to be computed in both intervals. Finally, the phase space portraits shown by the quantity consist of 622521 i.c. for the start space and 594690 i.c. for the start space.

In order to make appropriate comparisons on each start space, we need to evaluate the and the SElLCE on equal standings. Thus, it is necessary to decide a cut–off level for the distributions of their values to use the same scales for both indicators.

In Fig. 4 we show the histograms according to the final values of the techniques. The histograms on the top panel correspond to the start space and the histograms on the bottom panel correspond to the start space. A clear bimodal distribution (one peak corresponding to regular orbits, and the other to chaotic orbits) is shown for both indicators. In case of the start space, the indicators have less than values for the orbits below (pointed out on the top panel of the same figure). Thus, the scales on the start space go from to . The similarities between both histograms on the start space make easier the selection of a common cut-off level. However, there are important differences between the histograms on the start space. The SElLCE shows a bimodal distribution within the interval with its peaks well separated. The peaks in the distribution on the start space match the ones found for the same indicator, but on the start space (see both panels in Fig. 4). This is not the case for the . The peak that corresponds to the regular component is beyond the value (see Fig. 4, bottom panel and Fig. 2, bottom right panel). Nevertheless, we find that the number of orbits with values of the SElLCE below the value are negligible and the does not show further structure below that value. Thus, we decided to use that value as the cut–off level in the distributions on the start space to build the scales.

Notice also that, according to on the start space, the peak corresponding to regular orbits is . This value confirms the conclusions drawn in Section 4.1 with a sample composed only of 100 orbits for the same start space. That is, if we use the arithmetic means for the energy surface (Table 3), we arrive to a similar value of . Furthermore, the minimum between both peaks (which gives a clear distinction between the regular and the chaotic component) in the start space is taken around the value . This value also matches the one considered to separate regular from chaotic orbits for the same start space in Section 4.1.

In Fig. 5 we present the phase space portraits given by the SElLCE (left panels) and the (right panels) for both selected regions in the (top panels) and the (bottom panels) start spaces.

The phase space portraits of the region inside the start space given by the SElLCE and the present the same structures. In spite of this high level of coincidence in the description of the start space, there are differences in the start space. On the bottom right panel of Fig. 5, the quantity shows an amount of spurious structures10 not shown by the variational indicator (bottom left panel of the same figure). In fact, the spurious structures make difficult to select a threshold to distinguish between regular and chaotic orbits. Still, it is possible to define a cut-off value for the , but the procedure is more handmade. For example, if we take as regular orbits those preserving decimals in their fundamental frequencies (see Section 4.1 and the previous discussion), we recover the portrait given by the SElLCE for the start space. However, the SElLCE (as well as most of the VICs studied in previous papers, for instance M11 and D12) shows a large separation of the different kinds of motion in both portraits (see also Fig. 4), and the choice of a threshold (if it is not already defined) is easier.

These results, which are not as clear as those given by the variational indicators in terms of motion separation, make the FMFT a second choice to study the global dynamics in a divided phase space.

The procedure used to determine the chaoticity or regularity of the orbits with the FMFT is standard. The deficiency in the descriptions is basically due to the method’s sensitivity to its parameters. Therefore, the reliability of the FMFT as a global indicator of chaos is, comparatively speaking, limited.

#### The computing times

The VICs need to compute the equations of motion together with the variational equations to determine the nature of the orbits. This process may take a large amount of time. On the other hand, the FMFT works with the computation of the frequencies, which is a fast process, but the determination of the to characterise the nature of the orbits involves further calculations. For example, the computation of the SElLCE (a fast VIC) to obtain the top and bottom left panels of Fig. 5 took 670 hs and 330 hs, respectively, by using the DOPRI8 routine. Instead, the computation of the fundamental frequencies (nothing more than the the line of biggest amplitude in each d.o.f. was computed) with the FMFT for 624100 and 594690 i.c. took 30 hs for both experiments, after the integration of the equations of motion. The time spent in this integration for u.t. was of 570 hs and 280 hs for the regions in the and the start spaces, respectively. In both cases we used the Taylor method (see Section 3.4). Furthermore, if we want to use the FMFT as a global indicator of chaos, we need to compute the , i.e. we need to compute the frequency vector for two time intervals. If we overlap the time intervals in a 50%, we economize on computing time. The total time interval must be u.t. under these circumstances. That is, 50% larger than the one used by the variational indicators. Finally, the time required by the FMFT extends from 600 (570 + 30 hs) to 885 hs and from 310 (280 + 30 hs) to 450 hs to obtain the portraits on the top and bottom right panels of Fig. 5, respectively.

The advantage of the fast computation of the frequencies is certainly lost in the time–consuming process involved in the computation of the .

In conclusion, the FMFT as a global indicator of chaos seems to be rather inconvenient when there are other alternatives such as fast VICs like the SElLCE.

## 5 The selection of proper VICs

The aim of this section is to select the proper VICs to be used in the forthcoming experiments on the and the start spaces of the ScTS model (Section 6).

We first study a small region in the start space (within the intervals and ) on the energy surface , shown in the red square in Fig. 6). We use the MEGNO, the SElLCE, the FLI, the OFLI and the APLE. This first experiment provides us with a comparison between quite new techniques (the SElLCE and the APLE) and VICs with excellent performances shown in previous studies (the MEGNO and the FLI/OFLI, e.g. in M11 and D12).

In addition, we use the same characteristic times along the experiments.

We are interested in the resolving power and the speed of convergence of the VICs. Therefore, we compute 10201 i.c. within the region marked in Fig. 6 for u.t.

### 5.1 The resolving power of the VICs

First, we briefly discuss the resolving power of the VICs under the circumstances of the experiment.

As we are looking for VICs to study big samples of i.c., we consider the final values of the indicators rather than their time evolution curves. Furthermore, when the indicator grows or decreases exponentially, it is convenient to stop the integration process at a particular saturation value. Thus, instead of simply registering the final values of the indicator, we can record the time needed to reach such saturation value, i.e. the saturation times (Skokos et al., 2007). The final values together with the saturation times of the VIC should be considered for the description of the phase space portrait (M11 and D12).

In Fig. 7 we present the phase space portraits given by the final values of the MEGNO and the SElLCE (top left and top right panels, respectively); the OFLI final values (middle left panel) and the saturation times (middle right panel); and the APLE final values (bottom left panel) and the saturation times (bottom right panel). The phase space portraits given by the FLI are very similar to those shown with the OFLI and they are not included for the sake of simplicity.

The region under analysis has a regular component (the main central structure) and a chaotic component (surrounding the main central structure). Moreover, there are two different regions of chaos. The region that surrounds the regular component is composed of very strong chaotic orbits while the band on the right side of every panel in Fig. 7 is filled with moderate to strong chaotic orbits.

The values of the MEGNO and the SElLCE do not grow exponentially for chaotic orbits as it happens, for instance, with the OFLI (and the FLI) or the APLE. Then, a saturation value is not defined and a saturation time is not recorded in the computing process. However, the values of the MEGNO and the SElLCE can be very large for strong chaotic orbits and very different from their values for regular and mild chaotic orbits. In a divided phase space with regions of very strong chaos, the large values of these indicators for strong chaotic orbits may hide the small differences between regular and mild chaotic orbits. Indeed, the differences among strong chaotic orbits are generally not as important as the differences between strong and mild chaotic orbits or between mild chaotic orbits and regular orbits. Thus, a saturation value for the MEGNO or the SElLCE is also very useful in this kind of scenario. As shown in Fig. 6, we take the saturation value for the MEGNO and its associated value for the SElLCE11.

In case of the OFLI (and the FLI) or the APLE, the chaotic orbits grow exponentially fast and a saturation value is needed to avoid the propagation of errors in their computation. The saturation value in Fig. 7 for the OFLI (middle panels) is (see Section 3.2). The saturation value in Fig. 7 for the APLE (bottom panels) is (see Section 3.3).

In the top panels of Fig. 7, we can see the regular and the chaotic components well separated by the MEGNO (left panel) and the SElLCE (right panel). Furthermore, the two different chaotic regions are clearly shown: the region of strong chaos (surrounding the central structure with regular orbits) and the band on the right side of the panels filled with moderate to strong chaotic orbits. There is no structure at all in the main central structure composed of regular orbits (top left panel) because of the MEGNO’s asymptotic threshold (Section 3.1). Also, there is no further structure in the strong chaotic region due to the selected saturation value. The SElLCE shows a similar portrait (top right panel).

In middle panels of Fig. 7, we use the OFLI to describe the region under analysis. On the one hand, we recover the phase space portraits (middle left panel) given by the MEGNO and the SElLCE in the top panels of the same figure. On the other hand, the band on the right side of both panels shows a more detailed structure with the OFLI. This is due to the fact that the saturation value taken for the OFLI covers a wider range of chaotic orbits than the values chosen for the MEGNO or the SElLCE. The saturation times (middle right panel) enhance the information given with the final values of the OFLI, specially within the strong chaotic region.

In the bottom panels of Fig. 7, we present the phase space portraits given by the APLE final values (left panel) and saturation times (right panel). The description is very similar to those given by the other VICs if we consider both the final values and the saturation times of this indicator. For example, the identification of the band on the right side of the panels needs the information of the saturation times. The saturation times show a clear distinction between the region of strong chaos surrounding the main central structure and the region with mild and strong chaotic orbits in the band on the right side of the panels. This is not the case for the other VICs (the MEGNO, the SElLCE, the OFLI and the FLI) for which the saturation times play a complementary role and the final values were enough for such identification. Nevertheless, the degree of coincidence in the descriptions of the four VICs (plus the FLI) seems to be high. Thus, there is no decisive advantage in favour of one of them.

### 5.2 The speed of convergence of the VICs: thresholds

Now we focus on another important characteristic of the VICs which is fundamental for an efficient study of large samples of i.c.: the speed of convergence.

We use the same total integration time ( u.t.) and the same i.c. (10201) used in Section 5.1. The initial u.t. of the time interval does not give reliable information for our purposes and it is not considered in this experiment. The remaining u.t. are divided into sub-intervals of u.t. each. We compute the number of chaotic orbits in each sub-interval. Before that, we had to identify the chaotic orbits. So, we needed accurate thresholds for each VIC. Many of them have theoretical or empirical estimations of their thresholds. Yet, they are only approximations and should be adjusted to each particular situation. Hence, we generalize the experiments shown in M11 and D12 and consider a wide range of thresholds for the different VICs. Finally, for each sub-interval, we compute the number of chaotic orbits for every VIC according to the different thresholds.

Performing several tests (i.e. considering different thresholds and VICs), we find that a percentage of for the chaotic component is stable and prevails in the experiment. We take it as the “true” percentage for the chaotic component. Having such a “true” percentage, we can adjust the thresholds of the VICs under consideration. The fastest convergency of the VIC to a stable percentage of chaotic orbits close to the “true” percentage gives the most appropriate threshold for that particular indicator.

In order to identify such appropriate thresholds for every VIC under consideration, we define the “selected thresholds”. These thresholds allow the associated indicator to reach the “true” percentage of chaotic orbits within the total integration time12.

In Fig. 8 we present the variation of the selected thresholds with the integration time for the APLE (top left panel), the MEGNO (top right panel), the FLI (bottom left panel) and the OFLI (bottom right panel). The SElLCE is not included in the experiment because it has not a defined mechanism to compute its threshold. The filled circles in Fig. 8 indicate the selected thresholds that allow the convergency of the corresponding VIC to a stable, close to the “true” percentage of chaotic orbits. The filled squares in the same figure indicate the selected thresholds that offer the fastest convergency. In other words, the filled squares indicate the most appropriate thresholds for a given indicator.

The APLE has the theoretical time–independent threshold (see Section 3.3 and Lukes-Gerakopoulos et al. (2008) for further details). However, we consider threshold values between and . The threshold values below characterise 100% of the orbits as chaotic orbits for the whole time interval and thus, they are discarded. The threshold values above do not arrive at the “true” percentage of the chaotic component within the time interval and thus, they are also discarded. If we use the threshold values in the interval , the APLE converges to stable percentages of the chaotic component close to the “true” percentage (see top left panel of Fig. 8). If we use the threshold values around , the APLE’s convergency is the fastest (approximately by u.t.).

The MEGNO, like the APLE, has also a time–independent threshold: (see Section 3.1 and Cincotta et al. (2003) for further details). For the same reasons given for the APLE, we consider threshold values above the theoretical approximation, i.e. between and . The threshold values below characterise 100% of the orbits as chaotic orbits for the whole time interval and thus, they are discarded. The threshold values above are much larger than the theoretical asymptotic estimation of and thus, they are also discarded. If we use the threshold values close to , the MEGNO reaches percentages of the chaotic component close to the “true” percentage more rapidly (see top right panel of Fig. 8). Nevertheless, using this range of selected thresholds the percentages are not convergent. Moreover, the percentages oscillate within the interval depending on the threshold selected. In other words, the MEGNO gives a higher percentage of chaotic orbits than the prevailing percentage shown by the other VICs. This experiment confirms the results yielded in M11 and D12.

The FLI and the OFLI have both similar behaviour and time–dependent thresholds to distinguish chaotic motion from regular motion (see Section 3.2 and Fouchard et al. (2002) for further details). After computing the logarithm in their definitions, a tentative threshold for them is with the integration time. We consider the general expression in order to study a variation of the thresholds changing the values of the parameter .

The values of the parameter for the FLI belong to the interval . The theoretical threshold with is not efficient because it does not arrive at stable percentages of the chaotic component close to the “true” percentage for the whole time interval. The same applies to the thresholds with so they are discarded. Furthermore, the threshold values with do not reach at the “true” percentage of the chaotic component within the time interval and thus, they are also discarded. If we use the threshold values in the interval in the time interval u.t., the FLI converges to stable percentages of the chaotic component close to the “true” percentage (see the bottom left panel of Fig. 8). If we use the threshold values with , the FLI’s convergency is the fastest (approximately by u.t., i.e. in the figure). These results are very similar to those obtained earlier with the APLE. Yet, there is an important difference: the threshold selected for the FLI () is closer to the theoretical approximation () than the threshold selected for the APLE (). Therefore, this indicates that the FLI may have a better theoretical approximation of its threshold than the APLE.

Finally, we also observe similar results for the OFLI. However, the parameter varies from to , i.e., for the first time the interval of selected thresholds includes the theoretical estimation (). That is, if we use the theoretical approximation, the OFLI reaches a stable percentage of chaotic orbits close to the “true” percentage of . Nevertheless, the convergency arrives at the end of the time interval, i.e. u.t. The threshold values below characterise 100% of the orbits as chaotic orbits for the whole time interval and thus, they are discarded. The threshold values above do not arrive at the “true” percentage of the chaotic component within the time interval and thus, they are also discarded. If we use the thresholds with in the interval and within the time interval u.t., the OFLI converges to stable percentages of the chaotic component close to the “true” percentage (see the bottom right panel of Fig. 8). If we use the threshold values with , the OFLI’s convergency is the fastest (approximately by u.t., i.e. in the figure). These results are very similar to those obtained with the FLI as mentioned before. Yet, once again there is an important difference: the threshold selected for the OFLI () is closer to the theoretical approximation () than the threshold selected for the FLI (). Therefore, this indicates that the theoretical approximation of the thresholds of the OFLI may work better.

## 6 Global analysis of the start spaces

In this section we proceed to study some dynamical aspects of the ScTS model by means of the complementary use of the VICs selected in Section 5.

We consider samples of around one millon i.c. on the and the start spaces for four different energy surfaces. Even though we do not study the remaining five start spaces indicated in Schwarzschild (1993); Papaphilippou & Laskar (1998), we gain sufficient information in order to briefly discuss some main features of the dynamics of the system.

The total integration times extend to u.t. for the energy surface (Section 4.1). The CPU times become critical and will be considered a priority in the selection of the VICs. The FLI family of indicators presents the most versatile indicators, showing good performances in all the experiments (see also M11 and D12 in which they are included in the corresponding CIsF). The MEGNO is also another candidate because of its good performance (also included in the CIsF of M11 and D12). The SElLCE seems to be a good alternative to the MEGNO when studying big samples of orbits. Furthermore, the MEGNO and the SElLCE show a lower computational cost in the experiment than that of the FLI family of indicators13. Therefore, we select the LI (the least time–consuming indicator given a fixed total integration time, see e.g. D12) and the pair MEGNO/SElLCE as the VICs to study the ScTS model.

### 6.1 The px0−pz0 start space

We have already mentioned that the SElLCE has certain advantages when we study big samples of orbits (see Section 3.1). Nevertheless, it does not have any defined procedure to determine a threshold in order to distinguish between regular and chaotic motion. Indeed, this fact is a serious problem if we want to separate regular motion from weak chaotic motion. We need a threshold for the SElLCE to correctly describe the start space of the ScTS model through the final values of the indicator. For the determination of such a threshold we require again a “true” percentage of chaotic orbits for the sample under study. The MEGNO shows high sensitivity to changes in its threshold (Section 5.2, M11 and D12) and hence it is not a reliable indicator to calibrate other VICs. Besides, the LI has a theoretical threshold (, with the total integration time), which is a still valid estimation as a decent first approximation. We start with the corresponding theoretical estimations14 and calibrate them using inspections of different profiles of i.c. within the sample. Then, we find reliable thresholds for the LI for every energy surface. We compute the “true” percentages of chaotic orbits with the LI and calibrate the thresholds of the SElLCE to best fit these values. In Table 4, we present the results ordered by energy surface. Notice that the thresholds associated with the SElLCE are considerably smaller than the thresholds of the LI, because the estimation of the lLCE given by the slope of the MEGNO converges more rapidly to the lLCE than the standard procedure (see Section 3.1 and Cincotta et al. (2003)).

A comment regarding the efficient way to work with the LI and the SElLCE is in order before going forward. On one hand the LI is a reliable, but slow indicator so that it should be applied for longer integration times but on smaller samples in order the computing times to be moderate. Then, we can estimate the percentages of chaotic orbits. On the other hand, we have the SElLCE, which has a fast speed of convergence so it can be applied for shorter integration times but on the whole sample of orbits (once again, the computing times are reasonably short). Finally, we calibrate the thresholds of the SElLCE adjusting the percentages of chaotic orbits to those obtained with the LI. This is an efficient way to take advantage of the different characteristics (reliability and speed of convergence) of both VICs. Nevertheless, the priority in this experiment is to obtain reliable percentages of chaotic orbits rather than saving computing time in order to obtain the more accurate phase space portraits. This is the reason why we previously applied both the LI and the SElLCE on the same samples and for the same total integration times.

Now that we have the thresholds for the SElLCE, we apply the latter to a sample of 1000444 i.c. to study the start space for the following four energy surfaces15: , , and (see Fig. 9).

In Fig. 9 we observe that independently of the energy surface, the chaotic component is the dominant one on the phase space portraits (yellow and orange colours. See also column three and five in Table 4). The regular component (red, blue and black colours) increases for lower energy values. However, this increment is not considerable yet.

On the top left panel of Fig. 9, we observe a clear separation of both components. We show a fully connected chaotic component for values of and the regular component for values of , besides some hyperbolic structures.

These structures grow as we move to lower energy surfaces. They are resonances that start to overlap with each other and occupy the regular component. We can also distinguish a separation inside the chaotic component. The fully connected chaotic domain moves back to lower values of . The region is occupied by another chaotic domain characterised by a resonance overlapping regime and a lower Lyapunov number (top right panel of Fig. 9). On the energy surface (bottom left panel of Fig. 9), the fully connected chaotic domain has noticeably moved backwards and the resonances occupy extended regions of the start space. Furthermore, a major resonance is seen around the value and minor resonances are seen around the values and . On the bottom right panel of Fig. 9, for the energy surface , the portrait shows how the resonances occupy most of the regular component. Moreover, major resonances grow within the chaotic domains. The most evident resonance is that located at value , which has already been identified in outer energy surfaces.

In spite of the short discussion of the start space given above, the efficiency of variational indicators to describe the global characteristics of the phase space is evident. The VICs easily identify the regular and the chaotic domains and their mutual interplay, along with phase space mixing processes like resonance overlapping. However, the identification of the underlying resonant structure by the frequency vectors is not readily understood. We should study the dimensionality of the tori (e.g. with the GALI method, Skokos et al. (2007), D12) to identify resonant orbits, locate the periodic orbits (e.g. with the OFLI, Fouchard et al. (2002), D12) and analyse their stability (e.g. with the MEGNO, Cincotta et al. (2008)) to look for resonant orbit families. The process could be time–consuming and rather inefficient since the VICs are not the most suitable tools for the task. In fact, the SAMs are far better options (see Section 1 and references therein).

The dynamics of the start space can be briefly described as a dominant chaotic component based mostly on a resonance overlapping regime. This phase mixing process produces a more homogeneous chaotic component for the energy surface .

For further studies concerning the start space of the ScTS model refer to Muzzio et al. (2005) and Cincotta et al. (2008).

Now, we apply the same package of techniques to the start space, where the dynamics is totally different since it is dominated by the regular orbits.

### 6.2 The x0−z0 start space

Here, we implement the same procedure applied in Section 6.1 to study the start space. We obtain the phase space portraits with the SElLCE for four energy surfaces, calibrating the thresholds with the assistance of the LI, and we further discuss the global dynamics of the start space.

In order to determine the percentages of the chaotic component, we proceed like we did in Section 6.1. We present the results ordered by energy surface on Table 5.

In Fig. 10, we present the phase space portraits of the start space16 for the same energy surfaces considered in Section 6.1. Also the total integration times remain the same and the VICs used for the study guarantee reliable and stable pictures. The numbers of orbits in the samples are different and depend on the energy surface: 998938 i.c. for the energy surface , 998948 for , 999400 for and 998626 for .

There is a significant change in the dynamics of the start space with respect to the results for the start space. The corresponding phase space portraits presented in Fig. 10 show a dominant regular component (see Table 5, columns three and five) instead of the rather large chaotic component seen for the start space. The percentages of the smaller component have a substantial increment toward the innermost energy surface. Indeed, in this case, the expansion of the chaotic component is due to the growth of the different orbital families separatrices. This fact is clearly seen in Fig. 10. On the top left panel of Fig. 10, we can see a very weak separatrix that separates the short–axis tubes from the outer long–axis tubes. On the top right panel of the same figure, we can see how the former separatrix is connected to the others which separate the inner long–axis tubes and the box orbits (see Papaphilippou & Laskar (1998) for further details on the different orbital families and their locations in phase space). This becomes evident in the bottom panels of Fig. 10. Furthermore, a strong chaotic domain appears in the crossing of the main separatrices. Therefore, as the dominant component comprises regular orbits, the main contribution to the chaotic component emerges from the separatrices of resonances.

Finally, the dynamics of the start space can be briefly described as a dominant regular component. The major resonances that separate the main orbital families have separatrices which grow when the negative energies increase. Furthermore, the growth of the separatrices increases the percentage of chaotic orbits very fast.

## 7 Conclusions

With this report we conclude a series of investigations (reported in Maffione et al. (2011a), M11 and D12) toward an efficient choice of a minimal package of techniques to study a general Hamiltonian.

Here, we apply both VICs and a SAM (the FMFT, Sidlichovský & Nesvorný (1997)) to study a fairly realistic dynamical model of an elliptical galaxy (the ScTS model, Muzzio et al. (2005)).

In order to select the appropriate VICs for the experiments, we not only consider previous comparisons, but also extend them with new indicators: the SElLCE and the APLE (Lukes-Gerakopoulos et al., 2008). The SElLCE is an efficient estimation of the lLCE by means of the MEGNO (Cincotta & Simó, 2000) and works as its reliable alternative in case of studying big samples of orbits. The APLE has not shown advantages over the FLI/OFLI for the purposes of the experiment, but it behaves similarly.

The SElLCE seems to improve the performance of the MEGNO, in particular situations. Therefore, it should be considered in the CIsF presented in D12 as a suitable alternative.

On the other hand, we consider the performance of the FMFT as the representative SAM to compare with the VICs. The FMFT is an improvement of the FMA outlined by Laskar (Laskar, 1990), which is widely used by the scientific community.

The main advantage of the SAMs is a fast computation of the frequencies. However, the performance of the FMFT as a global indicator of chaos is not as efficient as the performances shown by the VICs in the experiments. The SAMs are designed to compute the frequencies, which are quantities strictly related to regular motion. Using such techniques as global indicators of chaos involves forcing the methods to do something for which they are not designed. The VICs are based on the concept of local exponential divergence. Hence, the detection of chaos is their main purpose. A natural consequence of this is that an efficient application of those techniques is based on their complementary implementation (see Section 1 for references). The SAMs are of use to describe the resonance web while the VICs to study the interplay of regular and chaotic domains (see Section 6).

Finally, the recommended CIsF for the analysis of a general Hamiltonian is composed of the MEGNO/SElLCE, the FLI/OFLI and the GALI: the MEGNO/SElLCE and the FLI/OFLI as global indicators of chaos to obtain the phase space portraits and display the interaction between regular and chaotic orbits; the GALI to analyse small samples within regions of complex dynamics and regions where the chaotic component is dominant (see Skokos et al. (2007) and D12). Furthermore, the CIsF can be used with a SAM (e.g. the FMFT) in order to characterise the resonance web with the frequency vectors.

## Acknowledgments

The authors want to thank Dr. J.C. Muzzio for his assistance along the investigation and the English department of FCAG–UNLP for improving the English of the reported research. Also, we are grateful to the anonymous referee because of the insightful comments made to improve the original manuscript. This work was supported with grants from the Consejo Nacional de Investigaciones Científicas y Técnicas de la República Argentina (CCT–La Plata) and the Universidad Nacional de La Plata.

### Footnotes

1. footnotemark:
2. footnotemark:
3. footnotemark:
4. pagerange: Chaos detection tools: application to a self–consistent triaxial modelChaos detection tools: application to a self–consistent triaxial model
5. pubyear:
6. CIsF: “Chaos Indicators Function”. It is a minimal and efficient package of variational indicators to study a general Hamiltonian.
7. Let us remark that the following definition of the OFLI is slightly different from the definition of the OFLI given in Fouchard et al. (2002) where the authors use only one i.d.v.
8. There is a downloadable version of the FMFT in the personal web page of David Nesvorný: http://www.boulder.swri.edu/davidn/
9. The time intervals are and .
10. Some of them, are due to the phenomenon of Moiré, which is inherent in discrete Fourier transform techniques, Barrio et al. (2009b).
11. We can estimate the saturation value () for the SElLCE from the expression: , with u.t. (see Section 3.1).
12. Actually, the selected thresholds allow the associated VIC to reach a percentage of chaotic orbits between 82-84% within the total integration time.
13. In this test on computing times, we consider the definition of the FLI given in Froeschlé & Lega (2000), where the authors compute the evolution of the FLI using one deviation vector.
14. The theoretical estimations: (u.t.) for the energy surface ;