# [

## Abstract

We investigate the regular or chaotic nature of orbits of stars moving in the meridional plane of an axially symmetric galactic model with a dense, massive spherical nucleus and a dark matter halo component. In particular, we study the influence of the fractional portion of the dark matter, by computing in each case the percentage of chaotic orbits, as well as the percentages of orbits of the main regular resonant families. In an attempt to distinguish between regular and chaotic motion, we use the Fast Lyapunov Indicator (FLI) method to extensive samples of orbits obtained by integrating numerically the equations of motion as well as the variational equations. Furthermore, a technique which is based mainly on the field of spectral dynamics that utilizes the Fourier transform of the time series of each coordinate is used for identifying the various families of regular orbits and also to recognize the secondary resonances that bifurcate from them. Two cases are studied in our work: (i) the case where we have a disk galaxy model (ii) the case where our model represents an elliptical galaxy. Comparison with early related work is also made.

Influence of Dark Matter in Galaxies]Unveiling the Influence of Dark Matter in Axially Symmetric Galaxies
Caranicolas and Zotos]Nicolaos D. Caranicolas and Euaggelos E. Zotos^{1}

\jidPASA
\doi10.1017/pasa.2013.27
\jyear2018

haos – galaxies: kinematics and dynamics – galaxies: structure

## 1 Introduction

From time immemorial, the movements of stars have attracted the interest of the first astronomers. Several decades have been past, since astronomers have turned their interest to the kinematics of stars in galaxies and in galaxy clusters in general. Undoubtedly, galaxies are the basic units of the Universe. An interesting task is to determine the amount of matter in galaxies and, consequently, the amount of matter of the Universe. In particular, there are two key events that are closely associated with this interesting endeavour. The first is the pioneer research of J. Oort who studied the motion of stars in the neighborhood of the Sun and realized that in fact the stars were pulled by forces stronger than those that would be caused only by the visible matter (Oort, 1924). On the other hand, the astronomer F. Zwicky discovered that the fast motion of galaxies in clusters cannot be justified only by the visible mass (Zwicky, 1933). In other words, clusters of galaxies should be scrapped unless the their mass was much greater. The existence of dark matter is the most acceptable scenario to resolve these concerns. Things were clarified when Rubin & Ford (1970) presented observational evidence supporting the existence of dark matter in galaxies.

Knowing the rate ratio between luminous and dark matter in galaxies (disk or ellipticals), as well as the relations connecting the fundamental quantities which characterize both components is of great importance, since it allow us to understand how galaxies born and evolve and also how dark matter influences these procedures. Today, it is widely accepted that dark matter is the dominant element in galaxies, taking into account that the vast majority of the total matter of the Universe, 80% according to today measurements, is indeed dark. Dark matter seems to interacts through gravity. Moreover, apart form gravitational, no other type of dark matter interactions has been observed. This strongly indicates that dark matter interactions should be very weak, probably much more weaker than the particle physics weak interactions. The basic proposed candidates corresponding to dark matter are: (i) neutrinos (Hot Dark Matter), (ii) Warm Dark Matter (WDM), (iii) Cold Dark Matter (CDM) and finally (iv) weak interactive massive particles (WIMPS).

A strong indication for the presence of dark matter in galaxies is derived from their flattened rotation curves at large radii. Using Kepler’s third law we have , where is the rotational velocity at a radius , is the gravitational constant, while is the total mass within radius . Conducting observations at large galactocentric distances, where no luminous galactic component was present astronomers found, instead of declining at the expected rate , which is true if only , the velocity curves appear to be flattened. However, , implies that . This strongly suggests that the mass of galaxies continues to grow, even when there is no luminous component to account for this increase. Moreover, the profiles of the velocity curves indicate that the distribution of light does not match the distribution of mass. In other words, in each galaxy the mass to light ratio increases with the radius . The circular velocity needed for the construction of the rotation curve is obtained, as usually, by measuring the 21 cm emission line from neutral hydrogen HI (Rubin et al., 1980; Bosma, 1981; Clemens, 1985; Persic & Salucci, 1995; Honma & Sofue, 1997).

The existence of dark matter in elliptical galaxies has been confirmed using observational data derived either from hot gas or their X-ray emission (e.g., Loewenstein & White, 1999; Humphrey et al., 2006; Johnson et al., 2009; Das et al., 2010) or even using strong lensing methods (e.g., Rusin & Kochanek, 2005; Gavazzi et al., 2007; Koopmans et al., 2009; Faure et al., 2011). On the other hand, for the disk galaxies there is a plethora of scenarios describing the formation and the evolution of disk galaxies in correlation with dark matter (e.g., Dalcanton et al., 1997; Firmani et al., 1997; Avila-Reese et al., 1998; van de Bosch, 1998; Avila-Reese & Firmani, 2000). A large number of measurements of luminous and dark matter have become available over the last years (e.g., Cappellari et al., 2006, 2012; Wegner et al., 2012). Using these data in combination with dynamically modeling we can reconstruct and therefore study the orbital structure of galaxies.

Today we know that galaxies, contain large amounts of dark matter and, therefore, a further study could provide important information about this invisible matter. A simple way to do this, is to construct models of galaxies containing dark matter. The study of the dynamical properties of these models will provide interesting and useful information, which combined with data from observation will aid significantly in finding a solution for the problem of dark matter. The reader can find interesting information in the field of fitting mass models to the kinematics of disk and elliptical galaxies in a series of papers (see, e.g. Barnes et al., 2004; Chaktabarty, 2007; de Blok et al., 2008; Gebhardt & Thomas, 2009). Another important point, that needs to be emphasized, is the difference between the distribution of dark matter in galaxies and clusters. Observational data (e.g., Sahni, 2004) suggest that dark matter increases as we move away from the centre to the outer parts of the galaxies. On the contrary in galaxy clusters the dark matter distribution decreases with increasing distance from the galactic centre.

Even today, dark matter is still an open and controversial issue in Astronomy. This is true, because the standard model of cosmology (SMoC) seems to be incompatible with a large amount of data derived from extragalactic observations and modified gravity theories. The reader can find more interesting and detailed information regarding this issue in Kroupa (2012).

Taking into account all the above, there is no doubt, that dark matter plays an important role in the dynamical behavior of galaxies. On this basis, it seems of particular interest to build an analytical dynamical model describing the motion of stars both in disk and elliptical galaxies containing dark matter and then study, how the presence and the amount of dark matter affects the regular or chaotic nature of orbits as well as the behavior of the different families of orbits.

The present article is organized as follows: In Section 2 we present in detail the structure and the properties of our gravitational galactic model. In Section 3 we describe the computational methods we used in order to determine the character of orbits. In the following Section, we investigate how the parameter corresponding to the fractional portion of the dark matter in galaxies influences the character of the orbits, in both disk and elliptical galaxy models. Our paper ends with Section 5, where the discussion and the conclusions of this research are presented.

## 2 Description of the Galactic Model

In order to study the dynamical properties of galaxies astronomers often construct galactic models. A galactic model is usually a mathematical expression giving the potential or the mass density of the galaxy, as a function of the distance. The reader can find interesting models, describing motion in galaxies, in Binney & Tremaine (2008). Potential density pairs for galaxies were also presented by Vogt & Letelier (1996). Over the years, many galactic models have been proposed in order to model the orbital properties in axially symmetric systems. A simple yet realistic axisymmetric logarithmic potential was introduced in Binney (1981) for the description of galactic haloes at which the mass density drops like (see also Evans, 1993). However, the most well-known model for cold dark matter (CDM) haloes is the flattened cuspy Navarro-Frenk-White (NFW) model (Navarro et al., 1996, 1997), where the density at large radii falls like . This model being self-consistent has a major advantage and that’s why it is mainly used for conducting -body simulations. Our model, on the other hand, is not self-consistent but simple and contrived, in order to give us the ability to study in detail the orbital behavior of the galactic system. Nevertheless, contrived models can provide an insight into more realistic stellar systems, which unfortunately are very difficult to be studied if we take into account all the astrophysical aspects. Due to the fast that our gravitational model consists of Plummer type potentials, at large galactocentric distances the mass density decreases following the law.

In the present work, we shall investigate how the presence of the dark matter influences the character of the orbits in the meridional plane of an axially symmetric galaxy model with a spherical nucleus and a dark matter halo component. We shall use the usual cylindrical coordinates , where is the axis of symmetry.

The total potential in our model consists of three components: the main galaxy potential , the central spherical component and the dark matter halo component . The first one is represented by the new mass model introduced in Caranicolas (2012) and is given by

(1) |

where

(2) |

Here is the gravitational constant, is the mass of the galaxy, is the fractional portion of the dark matter in the galaxy, while represents the core radius of the halo of the disk. The shape of the galaxy is controlled by the parameters and which correspond to the horizontal and vertical scale length of the galaxy respectively. Therefore, this potential allow us to describe a variety of galaxy types from a flat disk galaxy when to an oblate elliptical galaxy when .

For the description of the spherically symmetric nucleus we use a Plummer potential (e.g., Binney & Tremaine, 2008)

(3) |

where and are the mass and the scale length of the nucleus, respectively. This potential has been used successfully in the past in order to model and therefore interpret the effects of the central mass component in a galaxy (see, e.g. Hasan & Norman, 1990; Hasan et al., 1993; Zotos, 2012a). At this point, we must make clear that Eq. (3) is not intended to represent the potential of a black hole nor that of any other compact object, but just the potential of a dense and massive nucleus thus, we do not include relativistic effects. The dark matter halo is modelled by a similar spherically symmetric potential

(4) |

where in this case, and are the mass and the scale length of the halo, respectively. The spherical shape of the dark halo is simply an assumption, due to the fact that galactic halos may have a variety of shapes.

In this work, we use the well known system of galactic units, where the unit of length is 1 kpc, the unit of mass is and the unit of time is yr. The velocity unit is 10 km/s, the unit of angular momentum (per unit mass) is 10 km kpc s, while is equal to unity. Finally, the energy unit (per unit mass) is 100 kms. In these units, the values of the involved parameters are: , , , and . For the disk model we choose: and , while for the elliptical model we have set , and . The fractional portion of dark matter , on the other hand, is treated as a parameter and its value varies in the interval .

It is well known, that in disk galaxies the circular velocity in the galactic plane ,

(5) |

is a very important physical quantity. A plot of for both disk and elliptical galactic model when is presented in Fig. 1(a-b), as a black curve. Moreover, in the same plot the green curve is the contribution from the luminous matter, while the red line corresponds to the contribution form the dark matter. It is seen, that in both cases, the component of the rotational curve generated by the dark matter remains flat for large galactocentric distances, while on the other hand, the luminous matter component continues to decline with increasing distance form the galactic center. We also observe, the characteristic local minimum of the rotation curve due to the massive nucleus, which appears at small values of , when fitting observed data to a galactic model (e.g., Irrgang et al., 2013; Gómez et al., 2010).

At this point, we must clarify, that the mass density in our new galaxy model obtains negative values when the distance from the centre of the galaxy described by the model exceeds a minimum distance , which strongly depends on the parameter . Fig. 2(a-b) shows a plot of vs for the both disk and elliptical galaxy models. We see, that even when the first indication of negative density occurs only when kpc, that is almost at the theoretical boundaries of a real galaxy. Here, we must point out, that our gravitational potential is truncated ar kpc for both reasons: (i) otherwise the total mass of the galaxy modeled by this potential would be infinite, which is obviously not physical and (ii) to avoid the existence of negative values of density.

Taking into account that the total potential is axisymmetric, the -component of the angular momentum is conserved. With this restriction, orbits can be described by means of the effective potential

(6) |

The term represents a centrifugal barrier; only orbits with small are allowed to pass near the axis of symmetry. The three-dimensional (3D) motion is thus effectively reduced to a two-dimensional (2D) motion in the meridional plane , which rotates non-uniformly around the axis of symmetry according to

(7) |

where of course the dot indicates derivative with respect to time.

The equations of motion on the meridional plane are

(8) |

while the equations describing the evolution of a deviation vector which joins the corresponding phase space points of two initially nearby orbits, needed for the calculation of the standard chaos indicators (the FLI in our case) are given by the variational equations

(9) |

Consequently, the corresponding Hamiltonian to the effective potential given in Eq. (6) can be written as

(10) |

where and are momenta per unit mass, conjugate to and respectively, while is the numerical value of the Hamiltonian, which is conserved. Therefore, an orbit is restricted to the area in the meridional plane satisfying .

## 3 Computational Methods

In order to obtain the mass profiles of galaxies, we have to construct dynamical models describing the main properties of the galaxies. These models can be generated by deploying two main techniques: (i) using superposition of libraries of orbits (e.g., Gebhardt et al., 2003; Thomas et al., 2004; Cappellari et al., 2006) or (ii) using distributions functions (e.g., Dejonghe et al., 1996; Gerhard et al., 1998; Kronawitter et al., 2000). In the literature there are also other more specialized dynamical models combining kinematic and photometric data. For instance, axially symmetric Schwarzschild models were used by Bridges et al. (2006), while Hwang et al. (2008) used Jeans models in order to fit observational data in the X-ray potential introduced by Humphrey et al. (2006). Furthermore, axisymmetric Schwarzschild models were also used by Shen & Gebhardt (2010) to fit data derived from the Hubble Space Telescope (HST).

In the recent years, Schwarzschild’s superposition method (Schwarzschild, 1979) has been heavily utilized by several researchers (e.g., Rix et al., 1997; Gebhardt et al., 2003; Thomas et al., 2004; Valluri et al., 2004; Krajnović et al., 2005; Thomas et al., 2005) in order to study dark matter distributions in elliptical galaxies therefore, we deem it is necessary to recall and describe briefly in a few words the basic points of this interesting method. Initially, closed cells define the configuration space, while orbits extracted from a given mass distribution construct the phase space. Then, integrating numerically the equations of motion, we calculate the amount of time spent by each particular orbit in every cell. Thus, the mass of each cell is directly proportional to the total sum of the stay times of orbits in every cell. Using this technique we manage to compute the unknown weights of the orbits, assuming they are not negative.

In our study, we want to know whether an orbit is regular or chaotic. Several chaos indicators are available in the literature; we chose the FLI method. The FLI (Froeschlé et al., 1997; Lega & Froeschlé, 2001) has been proved a very fast, reliable and effective tool, which is defined as

(11) |

where is a deviation vector. For distinguishing between regular and chaotic motion, we need to compute the FLI for a relatively short time interval of numerical integration . In particular, we track simultaneously the time-evolution of the orbit itself as well as the deviation vector in order to compute the FLI. The variational equations (9), as usual, are used for the evolution of the deviation vector. The main advantage of the FLI method is that only one deviation vector is required to be computed, while in the case of other dynamical indicators like SALI (Skokos et al., 2004) or GALIs (Skokos et al., 2007) more than one deviation vectors are needed. Therefore, by using the FLI method we need considerable less computation time for integrating and classifying massive sets of initial conditions of orbits.

The particular time-evolution of the FLI allow us to distinguish fast and safely between regular and chaotic motion as follows: when an orbit is regular the FLI exhibits a linear increase, while on the other hand, in the case of chaotic orbits the FLI evolution is super-exponential. The time evolution of a regular (R) and a chaotic (C) orbit for a time period of time units is presented in Fig. 3. We observe, that both regular and chaotic orbits exhibit an increase in the value of FLI but with a complete different rate. Unfortunately, this qualitative criterion is applicable only when someone wants to check the character of individual orbits by plotting and then inspecting by eye the evolution of FLI. Nevertheless, we can easily overcome this drawback by establishing a numerical threshold value, in order to quantify the FLI. After conducting extensive numerical experiments in different types of dynamical systems, we conclude that a safe threshold value for the FLI taking into account the total integration time of time units is the value 10. The horizontal, blue, dashed line in Fig. 3 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 FLI has become greater than the established threshold value. Therefore, if FLI the orbit is chaotic, while if FLI the orbit is regular.

In order to investigate the orbital properties (chaoticity or regularity) of the dynamical system, we need to establish some sample of initial conditions of orbits. The best approach, undoubtedly, would be have been to extract these sample of orbits from the distribution function of the model. Unfortunately, this is not available so, we followed another course of action. For determining the chaoticity of our models, we chose, for each set of values of the parameters of the potential, a dense grid of initial conditions in the phase plane, regularly distributed in the area allowed by the value of the energy . The points of the grid were separated 0.1 units in and 0.5 units in direction. For each initial condition, we integrated the equations of motion (8) as well as the variational equations (9) with a double precision Bulirsch-Stoer algorithm (e.g., Press et al., 1992). In all cases, the energy integral (Eq. (10)) was conserved better than one part in , although for most orbits, it was better than one part in .

Each orbit was integrated numerically for a time interval of time units (10 billion yr), which corresponds to a time span of the order of hundreds of orbital periods but of the order of one Hubble time. The particular choice of the total integration time is an element of paramount importance, especially in the case of the so called “sticky orbits” (i.e., chaotic orbits that behave as regular ones during long periods of time). A sticky orbit could be easily misclassified as regular by any chaos indicator^{2}

A first step towards the understanding of the overall behavior of our galactic system is knowing the regular or chaotic nature of orbits. Of particular interest, however, is also the distribution of regular orbits into different families. Therefore, once the orbits have been characterized as regular or chaotic, we then further classified the regular orbits into different families, by using the frequency analysis of Carpintero & Aguilar (1998); Muzzio et al. (2005). Initially, Binney & Spergel (1982, 1984) proposed a technique, dubbed spectral dynamics, for this particular purpose. Later on, this method has been extended and improved by Carpintero & Aguilar (1998) and Šidlichovský & Nesvorný (1996). In a recent work, Zotos & Carpintero (2013) the algorithm was refined even further, so it can be used to classify orbits in the meridional plane. In general terms, this method calculates the Fourier transform of the coordinates of an orbit, identifies its peaks, extracts the corresponding frequencies and search for the fundamental frequencies and their possible resonances. Thus, we can easily identify the various families of regular orbits and also recognize the secondary resonances that bifurcate from them.

Before closing this Section, we would like to make a short note about the nomenclature of orbits. All the orbits of an axisymmetric potential are in fact three-dimensional (3D) loop orbits, i.e., orbits that rotate around the axis of symmetry always in the same direction. However, in dealing with the meridional plane the rotational motion is lost, so the path that the orbit follows onto this plane can take any shape, depending on the nature of the orbit. We will call an orbit according to its behaviour in the meridional plane. Thus, if for example an orbit is a rosette lying in the equatorial plane of the axisymmetric potential, it will be a linear orbit in the meridional plane, etc.

## 4 Orbit Classification

In this Section, we will numerically integrate several sets of orbits in an attempt to distinguish the regular or chaotic nature of motion. We use the initial conditions mentioned in Sec. 3 in order to construct the respective grids of initial conditions, taking always values inside the Zero Velocity Curve (ZVC) defined by

(12) |

In all cases, the value of the angular momentum of the orbits is . We chose for both disk and elliptical galaxy models such energy levels which correspond to kpc, where is the maximum possible value of on the phase plane. Once the values of the parameters were chosen, we computed a set of initial conditions as described in Sec. 3 and integrated the corresponding orbits calculating the value of FLI and then classifying the regular orbits into different families.

### 4.1 Disk galaxy model

Figure | Type | |||
---|---|---|---|---|

4a | box | 2.09000000 | 0.00000000 | - |

4b | 2:1 banana | 5.79009906 | 0.00000000 | 1.83668227 |

4c | 1:1 linear | 5.47170604 | 30.97448116 | 1.23526267 |

4d | 3:2 boxlet | 1.58031642 | 16.62648333 | 3.54124377 |

4e | 4:3 boxlet | 11.37966251 | 0.00000000 | 4.83447045 |

4f | chaotic | 0.40000000 | 0.00000000 | - |

For the disk galaxy models we choose the energy level which is kept constant. Our investigation reveals that in our disk galaxy model there are six main types of orbits: (a) box orbits, (b) 1:1 linear orbits, (c) 2:1 banana-type orbits, (d) 3:2 resonant orbits, (e) 4:3 resonant orbits and (f) chaotic orbits. Note that every resonance is expressed in such a way that is equal to the total number of islands of invariant curves produced in the phase plane by the corresponding orbit. In Fig. 4(a-f) we present an example of each of the five basic types of regular orbits, plus an example of a chaotic one. In all cases, we set . The orbits shown in Figs. 4a and 4f were computed until time units, while all the parent periodic orbits were computed until one period has completed. The black thick curve circumscribing each orbit is the limiting curve in the meridional plane defined as . Table 1 shows the types and the initial conditions for each of the depicted orbits; for the resonant cases, the initial conditions and the period correspond to the parent periodic orbits.

To study how the fractional portion of the dark matter influences the level of chaos, we let it vary while fixing all the other parameters of our disk galaxy model. As already said, we fixed the values of all the other parameters and integrate orbits in the meridional plane for the set . In all cases, the energy was set to and the angular momentum of the orbits was . Once the values of the parameters were chosen, we computed a set of initial conditions as described in Sec. 3 and integrated the corresponding orbits computing the FLI of the orbits and then classifying regular orbits into different families.

In Figs. 5(a-f) we present six grids of orbits that we have classified for different values of the fractional portion of the dark matter . Here, we can identify all the different regular families by the corresponding sets of islands which are formed in the phase plane. In particular, we see the five main families already mentioned: (i) 2:1 banana-type orbits surrounding the central periodic point; (ii) box orbits are situated mainly outside of the 2:1 resonant orbits; (iii) 1:1 open linear orbits form the double set of elongated islands in the outer parts of the phase plane; (iv) 3:2 resonant orbits form the double set of islands above the box orbits; and (v) 4:3 resonant orbits correspond to the outer triple set of islands shown in the phase plane. Apart from the regions of regular motion, we observe the presence of a unified chaotic sea which embraces all the islands of stability. The outermost black thick curve is the ZVC defined by Eq. (12).

Fig. 6 shows the resulting percentages of the chaotic orbits and of the main families of regular orbits as varies. It can be seen, that there is a strong correlation between the percentage of most orbits and the value of . As the portion of the dark matter increases, there is a gradual decrease in the percentage of chaotic orbit, although this tendency is reversed in models with significant amount of dark matter . In particular, we observe that always chaotic motion is the dominant type of motion and when the amount of chaotic orbits grows at the expense of box orbits and 1:1 linear orbits. The meridional 2:1 banana-type orbits, on the other hand, are almost unperturbed by the shifting of the portion of dark matter. Moreover, the 4:3 resonant orbits exhibit a constant increase, while the percentage of the 3:2 resonant orbits remains at very low values. From the diagram shown in Fig. 6, one may conclude that the fractional portion of the dark matter affects mostly the 1:1, 4:3 resonant orbits and chaotic orbits in disk galaxy models.

Of particular interest, is to investigate how the variation in the fractional portion of the dark matter influences the position of the periodic points of the different families of periodic orbits shown in the grids of Fig. 5. For this purpose, we use the theory of periodic orbits (Meyer & Hall, 1992) and the algorithm developed and applied in ?. In Fig. 7(a-d) we present the evolution of the starting position of the parent periodic orbits of the four basic families of resonant orbits. The evolution of the 2:1 and 4:3 families shown in Figs. 7a and Fig. 7b respectively, is two-dimensional since the starting position of both families lies on the axis. On the contrary, studying the evolution of the 1:1 and 3:2 families of periodic orbits is indeed a real challenge due to the peculiar nature of their starting position . In order to visualize the evolution of these families, we need three-dimensional plots such as those presented in Figs. 7c and 7d, taking into account the simultaneous relocation of and . The stability of the periodic orbits can be obtained from the elements of the monodromy matrix as follows:

(13) |

where Tr stands for the trace of the matrix, and is called the *stability index*. For each set of value of , we first located, by means of an iterative process, the position of the parent periodic orbits. Then, using these initial conditions we integrated the variational equations in order to obtain the matrix , with which we computed the index . Our numerical calculations indicate, that in the disk galaxy models, all the different families of periodic orbits remain stable throughout the entire range of the values of .

### 4.2 Elliptical galaxy model

Figure | Type | |||
---|---|---|---|---|

8a | box | 1.90000000 | 0.00000000 | - |

8b | 2:1 banana | 7.18420419 | 0.00000000 | 2.47524161 |

8c | 1:1 linear | 1.53837422 | 33.45716659 | 1.48665183 |

8d | 3:2 boxlet | 1.06837848 | 22.64363815 | 4.44068043 |

8e | 4:3 boxlet | 12.92588226 | 0.00000000 | 5.92689152 |

8f | 5:3 boxlet | 1.58425772 | 9.55073196 | 7.37479109 |

8g | 8:5 boxlet | 12.42410864 | 0.00000000 | 11.83080833 |

8h | chaotic | 0.40000000 | 0.00000000 | - |

In the case of the elliptical galaxy model, we choose the energy level which is kept constant. Our numerical investigation shows that in our elliptical galaxy model there are seven main types of orbits: (a) box orbits, (b) 1:1 linear orbits, (c) 2:1 banana-type orbits, (d) 3:2 resonant orbits, (e) 4:3 resonant orbits, (f) 5:3 resonant orbits, (g) 8:5 resonant orbits and (h) chaotic orbits. It is worth noticing, that the basic resonant families, that is the 2:1, 1:1, 3:2 and 4:3 are common in both disk and elliptical galaxy models. However, in the case of the elliptical galaxy additional secondary resonances (i.e. 5:3 and 8:5) appears. Again, every resonance is expressed in such a way that is equal to the total number of islands of invariant curves produced in the phase plane by the corresponding orbit. In Fig. 8(a-h) we present an example of each of the seven basic types of regular orbits, plus an example of a chaotic one. In all cases, we set . The orbits shown in Figs. 8a and 8f were computed until time units, while all the parent periodic orbits were computed until one period has completed. Table 2 shows the types and the initial conditions for each of the depicted orbits; for the resonant cases, the initial conditions and the period correspond to the parent periodic orbits.

In order to study how the fractional portion of the dark matter influences the level of chaos, we let it vary while fixing all the other parameters in our elliptical galaxy model. As already said, we fixed the values of all the other parameters and integrate orbits in the meridional plane for the set . In all cases the energy was set to and the angular momentum of the orbits was . Once the values of the parameters were chosen, we computed a set of initial conditions as described in Sec. 3 and integrated the corresponding orbits computing the FLI of the orbits and then classifying regular orbits into different families.

Six grids of initial conditions that we have classified for different values of the fractional portion of the dark matter are shown in Figs. 9(a-f). By inspecting these grids, we can identify all the different regular families by the corresponding sets of islands which are produced in the phase plane. In particular, we see the seven main families of regular orbits already mentioned: (i) 2:1 banana-type orbits correspond to the central periodic point; (ii) box orbits situated mainly outside of the 2:1 resonant orbits; (iii) 1:1 open linear orbits form the double set of elongated islands in the outer parts of the phase plane; (iv) 3:2 resonant orbits form the double set of islands; (v) 4:3 resonant orbits correspond to the outer triple set of islands shown in the phase plane; (vi) 5:3 resonant orbits forming the set of the three small islands inside the region of box orbits and (vii) 8:5 resonant orbits producing the set of five islands. It is evident, that the structure of the phase plane in the elliptical galaxy models differs greatly from that of the disk models. We observe, that when the amount of dark matter in the elliptical galaxy is low, almost the entire phase plane is covered by different types of regular orbits. On the other hand, chaotic motion appears, mainly at the outer parts of the phase plane, only when the galaxy possesses a significant amount of dark matter ().

In the following Fig. 10 we present the resulting percentages of the chaotic orbits and of also the main families of regular orbits as varies. It can be seen, that the motion of stars in elliptical galaxies, is almost entirely regular, being the box orbits the all-dominant type. The percentage of box orbits is however reduced as the portion of the dark matter is increased, although they always remain the most populated family. It is also seen, that the percentages of the 2:1 banana-type orbits exhibits a minor decrease with increasing . On the other hand, the chaotic orbits start to grow rapidly as soon as the galaxy contains significant amount of dark matter (). Moreover, all the other resonant families of orbits are immune to changes of the portion of the dark matter, since their percentages remain almost unperturbed and at very low values (less than 5%). From Fig. 10 one may conclude that dark matter in elliptical galaxies mostly affects box and chaotic orbits. In fact, a portion of box orbits turns into chaotic as the galaxy gains more dark matter.

We close this Section, by presenting how the variation in the fractional portion of the dark matter influences the position of the periodic points of the different families of periodic orbits shown in the grids of Fig. 9(a-f). We follow the same method used previously in the case of the disk galaxy. In Fig. 11(a-f) the evolution of the starting position of the parent periodic orbits of the six basic families of resonant orbits is given. Once more, the evolution of the 2:1, 4:3 and 8:5 families shown in Figs. 11(a-c), is two-dimensional since the starting position of these families lies on the axis. On the contrary, the evolution of the 1:1, 3:2 and 5:3 families are shown in the three-dimensional plots in Figs. 11(d-f) thus following simultaneous the relocation of and as the fractional portion of the dark matter varies. Our numerical calculations suggest, that all the computed resonant periodic orbits were found to be stable. Furthermore, we should point out that several families of periodic orbits in the case of the elliptical galaxy do not cover the entire range of the values of .

## 5 Discussion and Conclusions

There is no doubt, that the determination of the nature of dark matter is one of the most interesting and challenging open problems that scientists try to solve. In the present paper, we used the analytic, axisymmetric, mass model which was introduced in Caranicolas (2012) and embraces the general features of a dense, massive nucleus and a spherical dark matter halo. We made this choice because observations show that the assumption of a spherical halo seems to be reasonable. On the other hand, non spherical and triaxial haloes are also possible in some galaxies (see, e.g. Gentile et al., 2004; Trachternach et al., 2008; Vera-Ciro et al., 2011). A galaxy with a dark matter halo is undoubtedly a very complex entity and, therefore, we need to assume some necessary simplifications and assumptions in order to be able to study mathematically the orbital behavior of such a complicated stellar system. For this purpose, our model is simple and contrived, in order to give us the ability to study different aspects of the dynamical model. Nevertheless, contrived models can provide an insight into more realistic stellar systems, which unfortunately are very difficult to be studied, if we take into account all the astrophysical aspects. On the other hand, self-consistent models are mainly used when conducting -body simulations. However, this is entirely out of the scope of the present paper. Once again, we have to point out that the simplicity of our model is necessary; otherwise it would be extremely difficult, or even impossible, to apply the extensive and detailed dynamical study presented in this study. Similar gravitational models with the same limitations and assumptions were used successfully several times in the past in order to investigate the orbital structure in much more complicated galactic systems (see, e.g. Zotos, 2012b, 2013a).

In this work, we investigated how influential is the parameter corresponding to the portion of the dark matter on the level of chaos and also on the distribution of regular families among its orbits in both disk and elliptical galaxy models. The main results of our research can be summarized as follows:

(1). In disk galaxy models, the fractional portion of the dark matter affects mostly the 1:1, 4:3 resonant orbits and the chaotic orbits, while the effect on all the other resonant families is very weak compared to them. In particular, chaotic motion is always the prevailing type of motion but when the amount of dark matter is high enough, the amount of chaotic orbits grows at the expense of box orbits and 1:1 linear orbits.

(2). That portion of dark matter in elliptical galaxy models influences mainly the box and the chaotic orbits. In fact, box orbits are the dominant family when the amount of dark matter is low, but the percentage of chaotic orbits quickly grows as the dark matter is being accumulated, by collapsing the percentage of box orbits, although they always remain the most populated family.

(3). The percentage of the observed chaos in disk galaxy models with dark matter is significantly larger compared to that in elliptical galaxy models. This result coincide with similar conclusions obtained using different types of dynamical models in order to model disk and elliptical galaxies (see, e.g. Papadopoulos & Caranicolas, 2005; Zotos, 2011).

The outcomes of the present research are considered as a promising first step in the task of exploring the orbital structure in both disk and elliptical galaxy models containing dark matter. Taking into account that our results are encouraging, it is in our future plans to modify properly our dynamical model in order to expand our study in three dimensions.

## Acknowledgments

We would like to express our warmest thanks to the anonymous referee for the careful reading of the manuscript and for all the aptly suggestions and comments which allowed us to improve both the quality and the clarity of our paper.

### Footnotes

- thanks: Corresponding author’s E-mail: evzotos@physics.auth.gr
- Generally, dynamical methods are broadly split into two types: (i) those based on the evolution of sets of deviation vectors in order to characterize an orbit and (ii) those based on the frequencies of the orbits which extract information about the nature of motion only through the basic orbital elements without the use of deviation vectors.

### References

- Avila-Reese, V., Firmani, C., Hernández, X. 1998, ApJ, 505, 37
- Avila-Reese, V., Firmani, C. 2000, RevMexAA, 36, 23
- Barnes, E.I., Sellwood, J.A., Kosowsky, A. 2004, AJ, 128, 2724
- Binney, J. 1981, MNRAS, 196, 455
- Binney, J., Spergel, D. 1982, ApJ, 252, 308
- Binney, J., Spergel, D. 1984, MNRAS, 206, 159
- Binney, J., Tremaine, S. 2008, Galactic Dynamics, 2 ed., Princeton Univ. Press, Princeton, USA
- Bosma, A. 1981, AJ, 86, 1825
- Bridges, T., Gebhardt, K., Sharples, R., Faifer, F. R., Forte, J. C., Beasley, M. A., Zepf, S. E., Forbes, D. A., Hanes, D. A., Pierce, M. 2006, MNRAS, 373, 157
- Cappellari, M., Bacon, R., Bureau, M., Damen, M. C., Davies, R. L., de Zeeuw, P. T., Emsellem, E., FalcÃ³n-Barroso, J., KrajnoviÄ, D., Kuntschner, H., et al. 2006, MNRAS, 366, 1126
- Cappellari, M., McDermid, R. M., Alatalo, K., Blitz, L., Bois, M., Bournaud, F., Bureau, M., Crocker, A. F., Davies, R. L., Davis, T. A., et al. 2012, Nature, 484, 485
- Caranicolas, N.D. 2012, MNRAS, 423, 2668
- Carpintero, D.D., Aguilar, L.A. 1998, MNRAS, 298, 1
- Chakrabarty, D. 2007, MNRAS, 377, 30
- Clemens, D.P. 1985, ApJ, 295, 422
- Dalcanton, J.J., Spergel, D.N., Summers, F.J. 1997, ApJ, 482, 659
- Das, P., Gerhard, O., Churazov, E., Zhuravleva, I. 2010, MNRAS, 409, 1362
- de Blok, W.J.G., Walter, F., Brinks, E., Trachternach, C., Oh, S.H., Kennicutt, R.C. 2008, ApJ 136, 2648
- Dejonghe H., de Bruyne V., Vauterin P., Zeilinger W.W. 1996, A&A, 306, 363
- Evans, N.W. 1993, MNRAS, 260, 191
- Faure, C., Anguita, T., Alloin, D., Bundy, K., Finoguenov, A., Leauthaud, A., Knobel, C., Kneib, J.-P., Jullo, E., Ilbert, O., et al. 2011, A&A, 529, A72
- Firmani, C., Avila-Reese, V., Hérnandez, X. 1997, ASP Conf. Ser., 117, 424
- Froeschlé, Cl., Gonczi, R., Lega, E. 1997, Planet. Space Sci., 45, 881
- Gavazzi, R., Treu, T., Rhodes, J. D., Koopmans, L. V. E., Bolton, A. S., Burles, S., Massey, R. J., Moustakas, L. A. 2007, ApJ, 667, 176
- Gebhardt, K., Richstone, D., Tremaine, S., Lauer, T. R., Bender, R., Bower, G., Dressler, A., Faber, S. M., Filippenko, A. V., Green, R., et al. 2003, ApJ, 583, 92
- Gebhardt, D., Thomas, J. 2009, ApJ, 700, 1690
- Gentile, G., Salucci, P., Klein, U., Vergani, D., Kalberla P. 2004, MNRAS, 351, 903
- Gerhard O., Jeske G., Saglia R.P., Bender R. 1998, MNRAS, 295, 197
- Gómez, F., Helmi, A., Brown, A.G.A., Li, Y.-S. 2010, MNRAS, 408, 935
- Hasan H., Norman C.A. 1990, ApJ, 361, 69
- Hasan H., Pfenniger D., Norman C. 1993, ApJ 409, 91
- Honma, M., Sofue, Y. 1997, PASJ, 49, 453
- Humphrey P.J., Buote D.A., Gastaldello F., Zappacosta L., Bullock, J.S., Brighenti F., Mathews W.G., 2006, ApJ, 646, 899
- Hwang, H. S., Lee, M. G., Park, H. S., Kim, S. C., Park, J.-H., Sohn, Y.-J., Lee, S.-G., Rey, S.-C., Lee, Y.-W., Kim, H.-I. 2008, ApJ, 674, 869
- Humphrey P.J., Buote D.A., Gastaldello F., Zappacosta L., Bullock, J.S., Brighenti F., Mathews W.G. 2006, ApJ, 646, 899
- Irrgang, A., Wilcox, B., Tucker, E., Schiefelbein, L. 2013, A&A, 549, A137
- Johnson, R., Chakrabarty, D., O’Sullivan, E., Raychaudhury, S. 2009, ApJ, 706, 980
- Koopmans, L. V. E., Bolton, A., Treu, T., Czoske, O., Auger, M. W., BarnabÃ¨, M., Vegetti, S., Gavazzi, R., Moustakas, L. A., Burles, S. 2009, ApJ, 703, L51
- Krajnović, D., Cappellari, M., Emsellem, E., McDermid, R.M., de Zeeuw, P.T. 2005, MNRAS, 357, 1113
- Kronawitter A., Saglia R.P., Gerhard O., Bender R., 2000, A&A Supplementary Series, 144, 53
- Kroupa, P. 2012, PASA, 29, 395
- Lega, E., Froeschlé, Cl. 2001, CeMDA, 81, 129
- Loewenstein, M., White, R.E. 1999, ApJ, 518, 50
- Meyer, K.R., Hall, G.R. 1992, Introduction to Hamiltonian Dynamical Systems and the N-Body Problem, Springer-Verlag, New York
- Muzzio, J.C., Carpintero, D.D., Wachlin, F.C. 2005, CeMDA, 91, 173
- Navarro, J.F., Frenk, C.S., White, S.D.M. 1996, ApJ, 462, 563
- Navarro, J.F., Frenk, C.S., White, S.D.M. 1997, ApJ, 490, 493
- Oort, J.H. 1924, Note on the Difference in Velocity between Absolutely Bright and Faint Stars, Proceedings of the National Academy of Sciences, 10, 253
- Persic, M., Salucci, P. 1995, ApJS, 99, 501
- Papadopoulos, N.J., Caranicolas, N.D. 2005, Baltic Astronomy, 14, 253
- Press, H.P., Teukolsky, S.A, Vetterling, W.T., Flannery, B.P. 1992, Numerical Recipes in FORTRAN 77, 2nd Ed., Cambridge Univ. Press, Cambridge, USA
- Rix H., de Zeeuw P.T., Cretton N., van der Marel R.P., Carollo C.M., 1997, ApJ, 488, 702
- Rubin, V.C., Ford, W.K. 1970, ApJ, 159, 379
- Rubin, V.C., Ford, W.K., Thonnard, N. 1980, ApJ, 238, 471
- Rusin, D., Kochanek, C.S. 2005, ApJ, 623, 666
- Sahni, V. 2004, Lecture Notes in Physics, 653, 141
- Schwarzschild, M. 1979, ApJ, 232, 236
- Shen J., Gebhardt K. 2010, ApJ, 711, 484
- Šidlichovský, M., Nesvorný, D. 1996, CeMDA, 65, 137
- Skokos, Ch., Antonopoulos, Ch., Bountis, T.C., Vrahatis, M.N. 2004, J. Phys. A: Math. Gen., 37, 6269
- Skokos, Ch., Bountis, T.C., Antonopoulos, Ch. 2007, Physica D, 231, 30
- Thomas, J., Saglia, R.P., Bender, R., Thomas, D., Gebhardt, K., Magorrian, J., Richstone, D. 2004, MNRAS, 353, 391
- Thomas, J., Saglia, R.P., Bender, R., Thomas, D., Gebhardt, K., Magorrian, J., Corsini, E.M., Wegner, G. 2005, MNRAS, 360, 1355
- Trachternach, C., de Blok, W.J.G., Walter, F., Brinks, E., Kennicutt, Jr.R.C. 2008, AJ, 136, 2720
- van den Bosch, F.C. 1998, ApJ, 507, 601
- Valluri, M., Merritt, D., Emsellem, E. 2004, ApJ, 602, 66
- Vera-Ciro, C. A., Sales, L.V., Helmi, A., Frenk, C.S., Navarro, J.F., Springel, V., Vogelsberger, M., White, S.D.M. 2011, MNRAS, 416, 1377
- Vogt, D., Letelier, P. 2005, PASJ, 57, 871
- Wegner, G.A., Corsini, E.M., Thomas, J., Saglia, R.P., Bender, R., Pu, S.B. 2012, AJ, 144, 78
- Zotos, E.E. 2011, New Astronomy, 16, 391
- Zotos, E.E. 2012a, New Astronomy, 17, 576
- Zotos, E.E. 2012b, ApJ, 750, 56
- Zotos, E.E. 2013a, PASA, 30, 12
- Zotos, E.E. 2013b, Nonlinear Dynamics, 73, 931
- Zotos, E.E., Carpintero, D.D. 2013, CeMDA, 116, 417
- Zwicky, F. 1933, Helvetica Physica Acta, 6, 110