1 Introduction

Galaxies behave as collisionless systems. If ${\cal N}\gg 1$ is the number of stars in a galaxy, the collisionless condition allows to reduce the full ${\cal N}$-particle probability density to the single-particle density, $f(\vec{x},t)$, where $\vec{x}$ is the phase space vector of any star. The density f obeys the collisionless Boltzmann equation in which the interaction term is reduced to a smooth potential, $\nabla^2\phi(\vec{q};t)\equiv 4\pi G{\cal N}\int f(\vec{q},\vec{p};t){\rm d}^3p$, generated by the whole galaxy.

Of most importance are steady state solutions of Boltzmann's equation. To find out these solutions we need a priori some knowledge about the orbital structure supported by the static potential $\phi(\vec{q})$. As the strong Jeans theorem states, f should depend on the phase space variables only through the integrals of motion, I₁, I₂, I₃, associated to $\phi(\vec{q})$ (see Binney & Tremaine 1987, BT87 hereafter). The theorem makes sense when the field is integrable, that is, when the motion is completely regular (and quasiperiodic). In general this is not the case for realistic galactic potentials. If irregular motion occupies a comparatively large fraction of phase space then a different interpretation of the strong Jeans theorem is needed (see, for example, Merritt & Valluri 1996).

In the present work we consider a given potential $\phi(\vec{r})$, being the orbits in this field the building blocks of a galaxy. This is a rough approximation to, say, an elliptical galaxy since we are neglecting the presence of gas, tidal effects, rotation, the irregularities of the potential due to the irregular dynamics, etc. Within all possible choices of $\phi(\vec{r})$, one should take realistic models. This was done for instance by Merritt & Friedman (1996), Merritt & Valluri (1996), where they take a Dehnen's (1993) law for the density that, in turn, leads to a triaxial potential. Thus $\phi$ is expected to be a complicated function of the position or, as in the latter case, no analytical expression can be obtained for it. In this direction, to perform numerical investigations it is convenient, in a first systematic study, to simplify even more the problem by assuming the potential given by a simple and well-behaved function of the position. This will be our present approach.

Having the potential the efforts are devoted to obtain a picture of the global dynamics. That is, for any energy level to look for families of orbits related to the main periodic orbits and also to identify stochastic components. This is suitable for 2D systems because resonant tori mean periodic orbits but for 3D systems resonant tori mean, in general, 2D tori. For the stochastic regions it is relevant to have estimates of the time-scale for the manifestation of chaos in the orbital motion, the so-called Lyapunov time defined as the inverse of the largest Lyapunov characteristic number (LCN hereafter - see, for example, Reichl 1992). Stochastic components for which the Lyapunov time is much larger than the Hubble time can be regarded as regular for practical purposes.

To cope with the whole problem one could apply analytical tools developed for the study of Hamiltonian systems. But several questions restrict this approach. As an example let us mention that almost all the useful theorems deal with a special kind of systems, the so-called near-integrable Hamiltonian systems, and can be written in terms of action-angle variables. Its main part only depends on the actions, while the remainder that depends on both actions and angles, is considered as a small perturbation. This approximation is largely used for some problems in Celestial Mechanics but in Galactic Dynamics it seems that is not possible, in general, to perform the separation between an unperturbed Hamiltonian and a perturbation. Therefore this nice set of variables appears to be unsuitable to study, globally, the dynamics of galaxies. Yet, some attempts to compute such variables can be found in the literature (see for example McGill & Binney 1990; Binney & Spergel 1984; Papaphilippou & Laskar 1996 - PL96 hereafter). Hence in almost all cases the study should be done by means of numerical tools and the obtained results rest mainly in those numerical experiments rather than in rigorous theorems.

From the above discussion it turns out that it would be useful to have at hand a simple procedure to derive the basic dynamics. The phase space associated to a realistic Hamiltonian contains ordered and chaotic components. No additional global integrals besides the energy exist. However, if chaotic motion is confined to a comparatively small region of the phase space, regular motion should respect some other (perhaps local) "pseudo-invariants'' and then it is in principle possible to obtain a picture of the dynamics. In fact, this remark is what justifies the attempts to compute (numerically) the actions in certain galactic potentials. This was also the spirit behind the search for the so-called third integral of motion, the Hénon & Heiles (1964) famous paper being a pioneer work in this sense. Generically a "third integral'' (or convergent normal form or adelphic integral using an old fashioned name due to Whittaker 1917) does not exist but, for some ranges of the energy (or any other parameter), formal expansions can be useful as quasi-invariants to have a good description of the phase space. Different approaches to construct equilibrium models of galaxies were done by means of fully integrable potentials, like that of Stäckel form (see, for instance, Dejonghe & de Zeeuw 1988; de Zeeuw & Pfenniger 1988) or using models that allow for an approximate third integral (Petrou 1983a,b - see Sect. 2.1).

In recent years, however, the discussion is focused on the relevance of chaotic motion in models of elliptical galaxies. There is numerical evidence that, for instance, triaxial systems with a central strong cusp (simulating a mass concentration or black hole) contain a large amount of chaos (see Merritt & Valluri 1996; Valluri & Merritt 1998). A wide-spread tool used to identify regular and chaotic orbits and to estimate the time-scale over which chaos is relevant, is the computation of the LCN (or the KS entropy) for a set of orbits over motion times which are periods. Motion times about 10⁴ periods are lower bounds to reach a good estimate of the LCN. Indeed $\sim 10^5-10^6$ periods are necessary to obtain an accurate determination (see Sect. 4). For the sake of comparison and to keep in mind a realistic time scale, we recall that the Hubble time is in the order of 10², 10³ periods as much. In any case, motion times periods turn out to be too large (in a computational sense) when we deal with 10⁵ or 10⁶ orbits. The situation is even worst if the potential has not a simple analytical expression, the computing time being quite long in this case. For instance, Merritt & Valluri (1996) reported that, for their model, the time needed to integrate a single orbit and its variationals over 10⁴ periods to get the two positive Lyapunov exponents was about 80 minutes in a DEC Alpha 3000/700 workstation.

Alternative techniques were proposed to separate ordered and stochastic motion, to classify orbits in families, to describe the global structure of phase space, but not to get the LCN in shorter times. In Sect. 6 we shall resume this point together with some comparisons with the new technique here presented (MEGNO).

In Sect. 2 we discuss a simple, heuristic but effective, way to understand the different types of orbits in general 2D non-axisymmetric galactic potentials. Even though some qualitative ideas behind this approach were sketched, for instance, by Binney & Spergel (1982), Contopoulos & Seimenis (1990), Cincotta (1993), Cincotta et al. (1996), in this work we take advantage of the pendulum model to show that the main families of orbits, loop and box, respect, besides the energy, the same approximate invariant of motion. This approximation rests on the so-called ideal resonance problem largely discussed in the past (see for example Garfinkel 1966). The explicit derivation of this rough third integral leads to illustrative pictures of the orbital distribution at a given energy level being then, in principle, unnecessary the numerical calculation of the action integrals.

In Sect. 3 we sketch a representation of the global phase space useful to display the full dynamics, for a given energy level, on a single picture. It rests on the computation of a standard surface of section but instead of using a 2D plane we consider the natural manifold for the problem, which in this case corresponds to the 2D sphere. In Sect. 4 we introduce the Mean Exponential Growth factor of Nearby Orbits (MEGNO). This new tool has proven to be useful for studying global dynamics and succeeds in revealing the hyperbolic structure of phase-space, the source of chaotic motion. The MEGNO provides a measure of chaos that is proportional to the LCN, so that it allows to derive the actual LCN but in realistic physical times, $\sim 10^3$ periods. At least for the example discussed here, the MEGNO seems to provide more information about the dynamics than any other technique used before in this kind of applications.

As an example of non-axisymmetric galactic potential we study the 2D logarithmic model. We choose this potential to perform the numerical study because of three different reasons: i) a field with a simple analytic expression is easier to deal with; ii) it is, perhaps, a rough realistic model for the field acting on a star moving on the equatorial plane of a barlike galaxy or in the meridian plane of a spatial axisymmetric non-rotating elliptical galaxy; iii) it was largely studied, among others, by Binney & Spergel (1982, 1984), BT87, Miralda-Escudé & Schwarzschild (1989) - MS89 hereafter -, Lees & Schwarzschild (1992) and particularly by PL96 by means of the frequency map analysis. Therefore it allows to compare their results with the ones obtained here.

Finally, in Sect. 5, a survey is made of the most relevant features of the global dynamics for significative values of the energy and flatness of the potential.

A future paper will be devoted to the study of full 3D potentials by means of these or other simple tools.