1 INTRODUCTION
It is well known that galaxies are not observed as isolated entities, but often appear in binary, triple, or multiple stellar systems. On the contrary, just as there are clusters of stars, there are also clusters of galaxies. This leads to gravitational galactic interactions. The Milky Way–Magellanic Clouds system is a very indicative example of galactic interaction (see Gardiner, Sawa, & Fujimoto Reference Gardiner, Sawa and Fujimoto1994; Lin, Jones, & Klemola Reference Lin, Jones and Klemola1995; Weinberg Reference Weinberg1995, Reference Weinberg1998; Van der Marel Reference Van der Marel2001; Putman et al. Reference Putman, Staveley-Smith, Freeman, Gibson and Barnes2003; Connors et al. Reference Connors, Kawata, Maddison and Gibson2004; Bekki & Chiba Reference Bekki and Chiba2005). Another interesting stellar system of interacting galaxies is the Andromeda galaxy M31, with its two smaller companion galaxies M32 and NGC 205. The gravitational effects of the M32 galaxy on the spiral structure of M31 were studied by Byrd (Reference Byrd1978). This interesting triple stellar system was also investigated in an earlier work, using a self-consistent computer simulation code with revealing results (see Vozikis & Caranicolas Reference Vozikis and Caranicolas1994).
One of the main reasons for the transition from ordered to chaotic motion is the presence of massive objects in the central regions of the galaxies. Stars reaching the galactic center on highly eccentric redial orbits can be scattered out of the galactic plane displaying chaotic motion (see Sellwood & Moore Reference Sellwood and Moore1999). Furthermore, a central mass concentration can strongly perturb the stellar box orbits in elliptical galaxies which become chaotic (Merritt Reference Merritt1996). Observational indications suggest the presence of a strong central mass concentration by a very sharply rising rotation curve. A second reason for chaotic behaviour in the galactic motion is the presence of strong external perturbations (see Zotos Reference Zotos2011c). Increasing the magnitude of perturbations weakens or destroys the stability of orbits, thereby increasing the percentage of the stochasticity in the stellar system. Resonances are also responsible for chaotic motion (see Contopoulos & Grosbøl Reference Contopoulos and Grosbøl1986; Cincotta, Giordano, & Perez Reference Cincotta, Giordano and Perez2006). Of course, resonances are responsible not only for the presence of chaos in galactic stellar systems but also for the chaotic motion in the solar system (see Wisdom Reference Wisdom1987; Henrard & Caranicolas Reference Henrard and Caranicolas1990). The reader can find interesting information about the chaotic motion in galaxies and its connection with observational data in Grosbøl (Reference Grosbøl2002).
An extended and diffuse stellar halo envelops the Milky Way. Although only an extremely small fraction of the stars in the Solar neighborhood belong to this halo, they can be easily recognised by their extreme kinematics and metallicities. Stellar populations with these properties can now be followed to distances in excess of 100 kpc using luminous tracers such as RR Lyraes, blue horizontal branch stars, metal-poor giants, and globular clusters (e.g. Oort Reference Oort1926; Baade Reference Baade1944; Eggen, Lynden-Bell, & Sandage Reference Eggen, Lynden-Bell and Sandage1962; Searle & Zinn Reference Searle and Zinn1978; Vivas & Zinn Reference Vivas and Zinn2006; Morrison et al. Reference Morrison2009).
Low surface-brightness features seen in projection around other galaxies aid in the interpretation of the Milky Way's stellar halo and vice versa. Diffuse concentric ‘shells’ of stars on 100 kpc scales around otherwise regular elliptical galaxies have been attributed to accretion events (e.g. Schweizer Reference Schweizer1980; Quinn Reference Quinn1984). Recent surveys of M31 (e.g. Ferguson et al. Reference Ferguson, Irwin, Ibata, Lewis and Tanvir2002; Kalirai et al. Reference Kalirai2006; McConnachie et al. Reference McConnachie2009) have revealed an extensive halo (to ∼150 kpc) also displaying abundant substructure. The surroundings of other nearby Milky Way analogues are now being targeted by observations using resolved star counts to reach very low effective surface brightness limits, although as yet no systematic survey has been carried out to sufficient depth (e.g. Zibetti & Ferguson Reference Zibetti and Ferguson2004; McConnachie et al. Reference McConnachie2006; de Jong, Radburn-Smith, & Sick Reference de Jong, Radburn-Smith, Sick, Funes and Corsini2008; Barker et al. Reference Barker, Ferguson, Irwin, Arimoto and Jablonka2009; Ibata, Mouhcine, & Rejkuba Reference Ibata, Mouhcine and Rejkuba2009). A handful of deep observations beyond the Local Group suggest that stellar halos are ubiquitous and diverse (e.g. Sackett et al. Reference Sackett, Morrisoni, Harding and Boroson1994; Shang et al. Reference Shang1998; Malin & Hadley Reference Malin, Hadley, Merritt, Valluri and Sellwood1999; Martínez-Delgado et al. Reference Martínez-Delgado, Pohlen, Gabany, Majewski, Peñarrubia and Palma2009; Faúndez-Abans et al. Reference Faúndez-Abans, Reshetnikov, de Oliveira-Abans and Fernandes2009).
Stellar halos formed from the debris of disrupted satellites are a natural byproduct of hierarchical galaxy formation in the Λ cold dark matter (ΛCDM) cosmology. The entire assembly history of a galaxy may be encoded in the kinematics, metallicities, ages, and spatial distributions of its halo stars. Even though these stars constitute a very small fraction of the total stellar mass, the prospects are good for recovering such information from the halos of the Milky Way, M31, and even galaxies beyond the Local Group (e.g. Johnston, Hernquist, & Bolte Reference Johnston, Hernquist and Bolte1996; Helmi & White Reference Helmi and White1999). In this context, theoretical models can provide useful ‘blueprints’ for interpreting the great diversity of stellar halos and their various sub-components, and for relating these components to fundamental properties of galaxy formation models. Alongside idealised models of tidal disruption, ab initio stellar halo simulations in realistic cosmological settings are essential for direct comparison with observational data.
Therefore, it seems very challenging to construct a galactic gravitational model of three degrees of freedom (3D), in order to study and reveal the character of the orbits of the stars, in a binary stellar system composed of a primary quiet or active galaxy and a small satellite companion galaxy. The aim of the present paper is to study the character (regular or chaotic) of motion in the primary galaxy (hereafter the galaxy) under the perturbation caused by a satellite galaxy and to connect the degree of chaos with parameters such as the mass of the nucleus of the galaxy or the distance between the two galaxies. We shall also study the behaviour of the velocities near the central region and try to connect this behaviour with the scale length of the nucleus. Furthermore, we shall deal with the case where adiabatic mass transfer takes place in the galaxy, that is, the galaxy evolves with time. In this initial research, only the case where the two bodies move in circular orbits about the center of mass of the system is considered, for simplicity. We plan to study the case of the satellite in inclined, elliptic orbits in a future work.
In a recent paper (Zotos Reference Zotos2012b), we investigated in detail and revealed the regular or chaotic character of motion in a galactic gravitational model of three degrees of freedom describing a binary quasar system. According to the model, two quasars are hosted in a pair of interacting disk galaxies. Here we must emphasise that in the present research the second body of the binary stellar system, that is, the satellite companion galaxy, is considered and treated as a point mass. On the contrary, in Zotos (Reference Zotos2012b), we deal with a much more complicated binary system, since none of the two interacting disk galaxies can be treated as point mass due to the fact that both of them are hosting quasars in their cores. Thus, we may say that Zotos (Reference Zotos2012b) is a reasonable extension and generalization of the present paper. Taking into account that in both papers we study the properties of motion of the stars in similar, but in no case identical binary systems, it is natural and well expected to use the same philosophy regarding the basic setup of the gravitational models, the handling of the numerical calculations, and also the same standard methods of galactic dynamics and celestial mechanics.
The present article is organised as follows. In Section 2, we present our gravitational dynamical model, which describes the motion in a binary stellar system. In Section 3, we provide an analysis of the two degrees of freedom (2D) system, considering orbits in the galactic plane (z = 0). In Section 4, we study the character of motion in the 3D system, using different kinds of dynamical methods. In Section 5, we present an interesting theoretical analysis in order to support the numerically obtained results. In Section 6, we use a 3D time-dependent model in order to follow the evolution of orbits as mass is transported adiabatically from the disk to the nucleus of the primary galaxy. We close this paper with Section 7, where the conclusion is presented and a comparison between the present theoretical results with observational data is made. In the same section, we include statements concerning the astrophysical relevance of this research work.
2 PRESENTATION OF THE DYNAMICAL MODEL
Our gravitational model describes the dynamical properties of a binary stellar system which is composed of a primary galaxy and a small satellite companion galaxy. The primary galaxy (hereafter the galaxy) consists of a disk and a spherical nucleus and is described by the potential
\begin{eqnarray}
V_{\rm P}(r,z) &=& V_{\rm d}(r,z) + V_{\rm n}(r,z) \nonumber \\
&=& - \frac{G M_{\rm d}}{\sqrt{b^2 + r^2 + (a + \sqrt{h^2 + z^2})^2}} \nonumber \\
&&- \frac{G M_{\rm n}}{\sqrt{r^2 + z^2 + c_{\rm n}^2}}\,,
\end{eqnarray}
\begin{equation}
V_{\rm S}(r,z) = - \frac{G M_{\rm s}}{d},
\end{equation}
In our study, we shall use the theory of the circular restricted three-body problem. The two bodies move in circular orbits in an inertial frame OXYZ, with the origin at the center of mass of the system, with a constant angular frequency Ω p < 0, given by Kepler's third law
\begin{equation}
\Omega _p= \sqrt{\frac{G M_{\rm t}}{R^3}},
\end{equation}
\begin{equation}
\Phi _{\rm t}(x,y,z) = \Phi _{\rm P}(x,y,z) + \Phi _{\rm S}(x,y,z) + \Phi _{\rm rot}(x,y),
\end{equation}
\begin{eqnarray}
\Phi _{\rm P}(x,y,z) &=& - \frac{G M_{\rm d}}{\sqrt{b^2 + r_{a1}^2 + (a + \sqrt{h^2 + z^2})^2}} \nonumber \\
&&- \frac{G M_{\rm n}}{\sqrt{r_1^2 + c_{\rm n}^2}}, \nonumber\\[-7pt]
\\[-7pt]
\Phi _{\rm S}(x,y,z) &=& - \frac{G M_{\rm s}}{r_2}, \nonumber \\
\Phi _{\rm rot}(x,y) &=& -\frac{\Omega _p^2}{2}\left[\frac{M_{\rm s}}{M_{\rm t}}r_{a2}^2 + R_{\rm s} r_{a1}^2 - R^2 \frac{M_{\rm s}}{M_{\rm t}}R_{\rm s} \right], \nonumber
\end{eqnarray}
\begin{eqnarray}
r_{a1}^2 &=& (x - x_1)^2 + y^2, \nonumber \\
r_{a2}^2 &=& (x - x_2)^2 + y^2, \nonumber\\[-7pt]
\\[-7pt]
r_1^2 &=& r_{a1}^2 + z^2, \nonumber \\
r_2^2 &=& r_{a2}^2 + z^2, \nonumber
\end{eqnarray}
\begin{eqnarray}
x_1 &=& -\frac{M_{\rm s}}{M_{\rm t}}R, \nonumber \\
x_2 &=& R \left(1 - \frac{M_{\rm s}}{M_{\rm t}}\right) = R + x_1, \\
R_{\rm s} &=& 1 - \frac{M_{\rm s}}{M_{\rm t}}. \nonumber
\end{eqnarray}
The angular frequency Ω p is calculated as follows. The two bodies circulate around their common mass center of the system with angular frequencies Ω p1 and Ω p2 given by
\begin{eqnarray}
\Omega _{p1} &=& \sqrt{\frac{1}{x_1}\left(\frac{-{\rm d} V_{\rm S}(r)}{{\rm d}r}\right)_{r=R}}, \nonumber\\[-7pt]
\\[-7pt]
\Omega _{p2} &=& \sqrt{\frac{1}{x_2}\left(\frac{{\rm d} V_{\rm P}(r)}{{\rm d}r}\right)_{r=R}}. \nonumber
\end{eqnarray}
The above setup in the rotating frame has been applied successfully in several previous papers studying the nature of binary stellar systems (see Caranicolas & Papadopoulos Reference Caranicolas and Papadopoulos2009; Caranicolas & Zotos Reference Caranicolas and Zotos2009a; Caranicolas & Innanen Reference Caranicolas and Innanen2009; Zotos Reference Zotos2012b). Thus, we may conclude that the dynamical model described in Equations (4) and (5) is a good and also realistic model for a pair of interacting galaxies at the given mass ratio. Furthermore, in a recent paper (Zotos Reference Zotos2012b), we used the same treatment and setup, in order to study the dynamical properties in a binary system of two interacting galaxies. This could be considered a reasonable extension of the present study, where the second galaxy here is assumed to be only a satellite companion. The fact that the primary galaxy is sufficiently apart from its satellite allows us to assume that the tidal phenomena are very small and therefore negligible. In particular, our numerical calculations indicate that the tidal phenomena are significant enough when the distance between the centers of the two galaxies is R ≤ 0.95. In the present study, we assume that the range of this distance is 1.5 ≤ R ≤ 3 and, therefore, the tidal forces can be neglected.
In this rotating frame, the equations of motion are
\begin{eqnarray}
\ddot{x} &=& -\frac{\partial \Phi _{\rm t}}{\partial x} - 2 \Omega _p \dot{y}, \nonumber \\
\ddot{y} &=& -\frac{\partial \Phi _{\rm t}}{\partial y} + 2 \Omega _p \dot{x}, \\
\ddot{z} &=& -\frac{\partial \Phi _{\rm t}}{\partial z}, \nonumber
\end{eqnarray}
\begin{equation}
J=\frac{1}{2}\big(p_x^2 + p_y^2 + p_z^2 \big) + \Phi _{\rm t}(x,y,z) = E_J,
\end{equation}
All numerical outcomes of the present work are based on the numerical integration of the equations of motion (9), which was made using a Bulirsh–Stöer routine in fortran77, with double precision in all subroutines (see Numerical Recipes in FORTRAN, 2nd edn, in Press et al. Reference Press, Teukolsky, Vetterling and Flannery1992). The accuracy of the calculations was checked by the consistency of the Jacobi integral (10), which was conserved up to the 12th significant figure.
In this article, we shall use a system of galactic units where the unit of length is 20 kpc, the unit of mass is 1.8 × 1011 M⊙, and the unit of time is 0.99 × 108 yr. The velocity unit is 197 km s−1, while G is equal to unity (see Vozikis & Caranicolas Reference Vozikis and Caranicolas1992). In these galactic units, we use the values a = 0.15, b = 0.368, h = 0.00625, c n = 0.0125, and M s = 0.2. The values of the above quantities of the dynamical system remain constant during this research, while the values of M d, M n, and R are treated as parameters. The above numerical values of the constant dynamical quantities of the system secure positive density everywhere and free of singularities.
Figures 1(a–c) show the contours of the projections of the isopotential curves Φt(x, y, z) = EJ on the (x, y), (x, z) and (y, z) primary planes, respectively. The values of the parameters are M d = 2, M n = 0.08, R = 1.5, and Ω p = 0.821922. Lighter colors indicate higher values of EJ . Figure 2 shows the contours of the isopotential curves Φt(x, y, z) = EJ , on the (x, y) plane, with some additional details regarding the stability and the structure of the dynamical system. Moreover, L 1, L 2, L 3, L 4, and L 5 are the five Lagrange equilibrium points, while C 1 and C 2 are the centers of the two galaxies at a distance of R = 1.5. At these equilibrium points, we have
\begin{equation}
\frac{\partial \Phi _{\rm t}(x,y,z)}{\partial x} = \frac{\partial \Phi _{\rm t}(x,y,z)}{\partial y} = \frac{\partial \Phi _{\rm t}(x,y,z)}{\partial z} = 0.
\end{equation}
(a–c) Contours of the projections of the isopotential curves Φt(x, y, z) = EJ on the (x, y), (x, z), and (y, z) planes.
Contours of the projections of the isopotential curves Φt(x, y, z) = EJ on the (x, y) plane. The five Lagrange equilibrium points are indicated as L 1, L 2, L 3, L 4, and L 5, while C 1 and C 2 are the centers of the two galaxies at a distance of R = 1.5.
The rotation curve of the primary galaxy when M d = 2 and M n = 0.08 is shown as the black line in Figure 3. In the same plot, the red line is the contribution from the spherical nucleus, while the blue line is the contribution from the disk–halo component. We observe that at small distances from the primary galactic center r ≤ 0.1 length units, which equals 2 kpc, dominates the contribution from the spherical nucleus, while at larger distances r < 2 kpc, the disk–halo contribution is the dominant factor.
The rotation curve of the primary galaxy is shown as the black line. The red line is the contribution from the spherical nucleus, while the blue line is the contribution from the disk–halo potential. (Colors available only in the electronic version of the article).
It would be very illuminating to compute the mass density distribution ρ(x, y, z) of our binary system, derived from the total potential V T(r, z), using the Poisson's equation
\begin{eqnarray}
\rho (x,y,z) &=& \frac{1}{4 \pi } \nabla ^2 V_{\rm T}(x,y,z) \nonumber \\
&=& \frac{1}{4 \pi }\left(\frac{\partial ^2 V_{\rm T}}{\partial x^2} +\frac{\partial ^2 V_{\rm T}}{\partial y^2} + \frac{\partial ^2 V_{\rm T}}{\partial z^2} \right).
\end{eqnarray}
(a–c) Contours of the projections of the isodensity curves ρ(x, y, z) = const on the (x, y), (x, z), and (y, z) planes.
(a–c) 3D plots of the distribution of the total mass density ρ(x, y, z) on the (x, y), (x, z), and (y, z) planes.
(a–c) The evolution of the total mass density along the x, y, and z axes.
3 STRUCTURE OF THE 2D HAMILTONIAN SYSTEM
In this section, we shall investigate the properties of motion in the Hamiltonian system of two degrees of freedom. This can be derived from Equation (10) if we set z = pz = 0. Then, the corresponding Hamiltonian is
\begin{equation}
J_2 = \frac{1}{2}\big(p_x^2 + p_y^2 \big) + \Phi _{\rm t}(x,y) = E_{J2},
\end{equation}
Figures 7(a–b) show the (x, px ), y = 0, py < 0 Poincaré phase plane, for the motion of a star in the primary galaxy, which obtained by numerical integration of the equations of motion (9). In this case, the distance between the centers of the two galaxies is R = 1.5 and Ω p = 0.821922. The value of the Jacobi integral is E J2 = −2.6. This particular value of the Jacobi integral will remain constant for all our numerical experiments. In Figure 7(a) we have M d = 2 and M n = 0.08, that is, the primary galaxy has a dense and massive nucleus. We observe a large chaotic area, while the regular region is mainly confined around the stable retrograde periodic point. There are also some small islands of invariant curves embedded in the chaotic sea, corresponding to secondary resonant orbits. Note that near the center of the primary galaxy, the velocity obtains high values of the order of 900 km s−1. This results from the presence of the dense nucleus in the core of the galaxy and it is characteristic of strong galactic activity. The outermost black solid line defines the zero velocity curve. Figure 7(b) is similar to Figure 7(a), but when M d = 2.08 and M n = 0. As the total mass of the system is conserved, this means that now the disk occupies all the available mass of the primary galaxy. In this case, the pattern has two main differences from the pattern shown in Figure 7(a). The first difference is that the chaotic region is smaller and the second is that the velocity near the galactic core is smaller, about 470 km s−1, than that of Figure 7(a). Moreover, some secondary resonances appear in the area of the retrograde orbits. Some small sticky regions are also present in both cases (active or quiet galaxy). Therefore, one can say that the structure of chaos in galaxies is a result not only of galactic interaction but also of the nuclear galactic activity. In other words, the presence of the dense and massive nucleus in the galactic core of the primary galaxy affects drastically the character of motion of the stars.
The (x, px ) Poincaré phase plane when R = 1.5 and Ω p = 0.821922. (a) M d = 2 and M n = 0.08; (b) M d = 2.08 and M n = 0.
Figure 8(a) is similar to Figure 7(a) but when the distance between the centers of two the galaxies is R = 2 and Ω p = 0.533854. The majority of the orbits around the direct and retrograde periodic point are ordered. There is a significant chaotic area mainly near the central region of the galaxy. A high-velocity value is once again observed near the galactic core of the primary galaxy, while some secondary resonances are also present. Figure 8(b) is similar to Figure 8(a), but when the galaxy is quiet, that is, when M n = 0. Here, the secondary resonances look more prominent, while only a small and confined chaotic layer around the direct periodic point is present. Furthermore, as the primary galaxy is quiet, the velocity near the center of the galactic core is smaller than that observed in Figure 8(a). Figures 9(a–b) are similar to Figures 8(a–b), but when the distance between the centers of the two galaxies is R = 2.5 and Ω p = 0.381995. In Figure 9(a), we observe that there is a relatively small chaotic layer confined mainly near the central region of the galaxy, while the rest of the phase plane is covered by invariant curves corresponding to regular orbits circulating around the stable direct and retrograde periodic point. Figure 9(b) is similar to Figure 9(a) when the galaxy is quiet. All orbits seem to be regular. Chaotic motion was not observed and, if present, is negligible. Figures 10(a–b) are similar to Figures 9(a–b) but when the distance between the centers of the two galaxies is R = 3 and Ω p = 0.290593. In Figure 10(a), we observe that the entire phase plane is covered by invariant curves corresponding to regular orbits circulating around the stable direct and retrograde periodic point, while there is no indication of chaotic motion. Figure 10(b) is similar to Figure 10(a), but when the galaxy is quiet. Once again, all orbits are regular and chaotic motion was not observed. Therefore, we may conclude that the primary galaxy, for large separations from its companion, does not show chaotic motion irrespective of the presence of the dense and massive nucleus in the galactic core. The outermost solid curve in the (x, px ) phase planes shown in Figures 7–10 is the limiting curve in each case.
(a–b) Similar to Figures 7(a–b) but when R = 2 and Ω p = 0.533854.
(a–b) Similar to Figures 7(a–b) but when R = 2.5 and Ω p = 0.381995.
(a–b) Similar to Figures 7(a–b) but when R = 3 and Ω p = 0.290593.
Figure 11 shows the percentage A% on the (x, px ) phase plane covered by chaotic orbits as a function of the distance R between the centers of the two galaxies, in two different cases. In the first case, the primary galaxy is active and the values of the parameters are M d = 2 and M n = 0.08, while in the second case the primary galaxy is quiet and, therefore, the values of the parameters are M d = 2.08 and M n = 0. The value of the angular frequency Ω p is calculated for each particular value of the distance R from Equation (3). The range of the values regarding the distance between the centers of the two galaxies is 1.5 ≤ R ≤ 3. We observe in Figure 11 that the chaotic percentage A% decreases exponentially, as the distance R increases in both cases (active and quiet primary galaxy). Moreover, both fitting curves converge to zero, when R = 3. The lower fitting curve, which corresponds to the case when the primary galaxy is quiet, presents a more rapid reduction and reaches more quickly (when R ≥ 2.35) to the zero limit. A more detailed view of Figure 11 reveals that for each particular value of the distance R, the chaotic percentage A% is always smaller when the primary galaxy is quiet. We must point out that the chaotic percentage A% is calculated as follows: we choose 103 orbits with random initial conditions (x 0, p x0) in each phase plane and then divide the number of those who correspond to chaotic orbits to the total number of the tested orbits.
A plot of the area A% of the (x, p x ) phase plane covered by chaotic orbits as a function of the distance R, when the primary galaxy is active or quiet.
To have a better estimation about the degree of chaos in our dynamical system in every case, we have computed the average value of the Lyapunov characteristic exponent (LCE; see Lichtenberg & Lieberman Reference Lichtenberg and Lieberman1992). Figure 12 shows a plot of the 〈LCE〉 as a function of the distance R between the centers of the two galaxies, in two different cases. As previously, in the first case, the primary galaxy is active and the values of the parameters are M d = 2 and M n = 0.08, while in the second case the primary galaxy is quiet and, therefore, the values of the parameters are M d = 2.08 and M n = 0. The value of the angular frequency Ω p is calculated for each particular value of the distance R from Equation (3). The range of values regarding the distance between the centers of the two galaxies is again 1.5 ≤ R ≤ 3. One can see in Figure 12 that 〈LCE〉 decreases almost linearly as the distance R increases when the primary galaxy is active. On the contrary, when the nucleus is absent, we observe that the corresponding fitting curve presents a more rapid reduction. Once again, as pointed out in Figure 11, for each particular value of the distance R, the average value 〈LCE〉 is always smaller when the primary galaxy is quiet. Moreover, both fitting curves converge to zero when R = 3. Here we must note that it is well known that the value of the LCE is different in each chaotic component (see Saito & Ichimura Reference Saito, Ichimura, Casati and Ford1979). As we have in all cases regular regions and only one unified chaotic area in each (x, px ) phase plane, we calculated the average value of the LCE by taking 103 orbits with random initial conditions (x 0, p x0) in the chaotic area in each case and we integrated each orbit for a time interval of 105 time units in order to obtain a reliable value of the LCE regarding the nature of the orbit. Here we have to point out that all calculated LCEs corresponding to chaotic orbits were different only on the fifth decimal point in the same chaotic region.
A plot of the average value of the 〈LCE〉 as a function of the distance R, when the primary galaxy is active or quiet, for the 2D dynamical system.
In Figures 13(a–h) we present eight representative orbits of the 2D dynamical system. Figure 13(a) shows a quasi-periodic orbit circulating around the primary galaxy, with initial conditions x 0 = 0.167, y 0 = 0, and p x0 = 0, while the value of p y0 is obtained from the Jacobi integral (13) for all orbits. The values of all the other parameters are as in Figure 7(a). Figure 13(b) shows a periodic orbit moving around the center of the primary galaxy, with initial conditions x 0 = −0.1895, y 0 = 0, and p x0 = 0, while the values of all other parameters are as in Figure 7(b). This orbit is characteristic of the 1:2 resonance. In Figure 13(c), a figure-eight type periodic orbit, characteristic of the 1:3 resonance, is shown. The initial conditions are x 0 = 0.5788, y 0 = 0, and p x0 = 0, while the value of all other parameters are as in Figure 8(a). Figure 13(d) shows a periodic orbit moving around the primary galaxy with initial conditions x 0 = 0.4295, y 0 = 0, and p x0 = 0, while the values of all the other parameters are as in Figure 7(b). This orbit is a characteristic example of the 3:3 resonance. In Figure 13(e), a periodic orbit with initial conditions x 0 = 0.7445, y 0 = 0, and p x0 = 0 is presented. This orbit belongs to the family of the 2:3 resonant orbits. The values of all other parameters are as in Figure 7(a). In Figure 13(f), we see a periodic orbit with initial conditions x 0 = −0.9854, y 0 = 0, and p x0 = 0. The values of all other parameters for this orbit are as in Figure 7(b). Figure 13(g) shows a complicated resonant periodic orbit of higher multiplicity, circulating around the primary galaxy, with initial conditions x 0 = 0.673, y 0 = 0, and p x0 = 0, while the values of all the other parameters are as in Figure 7(a). This orbit produces a set of nine tiny islands of invariant curves which are embedded in the retrograde area of the unified chaotic domain shown in Figure 7(a). A chaotic orbit with initial conditions x 0 = −0.31, y 0 = 0, and p x0 = 0 is given in Figure 13(h). The values of all other parameters for this orbit are as in Figure 7(a). We observe that stars moving in chaotic orbits pass arbitrarily around the galactic core without getting close enough. All orbits shown in Figures 13(a–h) were calculated for a time period of 150 time units.
(a–h) Eight representative regular orbits of the 2D dynamical system. The values of the initial conditions and all the other parameters are given in the text.
4 STRUCTURE OF THE 3D HAMILTONIAN SYSTEM
In this section, we shall investigate the regular or chaotic nature of motion in the 3D Hamiltonian system described by Equation (10). In order to keep things simple, we shall use our experience gained from the study of the 2D dynamical system in order to obtain a clear picture regarding the properties of motion in the 3D dynamical model. We are particularly interested in locating the initial conditions in the 3D dynamical system, producing regular or chaotic orbits. A convenient way to obtain this is to start from the (x, px ) phase planes of the 2D system with the same value of the Jacobi integral as used in the 2D system and described in the previous section. Specifically, the regular or chaotic nature of the 3D orbits is found as follows: we choose initial conditions (x 0, p x0, z 0), y 0 = p z0 = 0, such that (x 0, p x0) is a point on the phase planes of the 2D system. The points (x 0, p x0) lie inside the limiting curve
\begin{equation}
\frac{1}{2}p_x^2 + \Phi _{\rm t}(x) = E_{J2},
\end{equation}
In Figure 14, we present a plot of the average value of the LCE of the 3D system, as a function of the distance R between the centers of the two galaxies, in two different cases. As in the previous section, in the first case, the primary galaxy is active and the values of the parameters are M d = 2 and M n = 0.08, while in the second case the primary galaxy is quiet and, therefore, the values of the parameters are M d = 2.08 and M n = 0. The value of the angular frequency Ω p is calculated for each particular value of the distance R from Equation (3). The range of values regarding the distance between the centers of the two galaxies is again 1.5 ≤ R ≤ 3. One can see in Figure 14 that the 〈LCE〉 decreases almost linearly as the distance R increases when the primary galaxy is active when R ≤ 2.25. On the contrary, when the nucleus is absent, we observe that the corresponding fitting curve presents a more rapid reduction. Once again, as pointed out in Figure 12 regarding the 2D system, for each particular value of the distance R, the average value 〈LCE〉 is always smaller when the primary galaxy is quiet. Moreover, both fitting curves converge earlier to zero value, than in Figure 12. The method we used in order to compute the 〈LCE〉 in each case is the same as described in Figure 12. For all tested chaotic orbits in the 3D system, the initial value of z 0 is common and equal to 0.1. Here we must point out that if we compare the plot of the 3D system shown in Figure 14 with that of the 2D system shown in Figure 12, we observe that in each case (active and quiet galaxy, respectively) the average values of the LCEs of the 3D system are smaller than those of the 2D system. Therefore, we may conclude that the degree of chaos in the primary galaxy is smaller when the dynamical system has three degrees of freedom.
A plot of the average value of the maximum 〈LCE〉 as a function of the distance R, when the primary galaxy is active or quiet, for the 3D dynamical system.
Our next step is to study the character of orbits with initial conditions (x
0, p
x0, z
0), y
0 = p
z0 = 0, such that (x
0, p
x0) is a point in the regular regions of Figures 7(a–b), 8(a–b), 9(a–b), and 10(a–b). The phase space of a conservative system of three degrees of freedom has six dimensions, i.e. in Cartesian coordinates
$\left(x, y, z, \dot{x}, \dot{y}, \dot{z}\right)$
. For a given value of the Jacobi integral, a trajectory lies on a five-dimensional manifold. In this manifold, the surface of the section is four-dimensional. This does not allow us to visualise and interpret directly the structure and the properties of the phase space in dynamical systems of three degrees of freedom. One way to overcome this problem is to project the surface of the section to space with lower dimensions. In fact, here, we will apply the method introduced by Pfenniger (Reference Pfenniger1984) (see also Revaz & Pfenniger Reference Revaz and Pfenniger2001). We take sections in the plane y = 0, py
< 0 of 3D orbits, whose initial conditions differ from the plane parent periodic orbits only by the z component. The set of the resulting four-dimensional points in (x, px
, z, pz
) phase space is projected on the (z, pz
) plane. If the projected points lie on a well-defined curve, we call it an ‘invariant curve’, then the motion is regular, while if not, the motion is chaotic. The projected points on the (z, pz
) plane show nearly invariant curves around the periodic points at z = 0, pz
= 0, as long as the coupling is weak. When the coupling is stronger, the corresponding projections in the (z, pz
) plane displays an increasing departure of the plane periodic point, up to making a direct orbit a retrograde one and vice versa. Here, we must define what one means by direct and retrograde 3D orbit. If consequents in the (z, pz
) section of the 3D orbit drop in one of the two domains of the corresponding section of 2D orbits at the same value of the Jacobi integral EJ
, we can distinguish between direct and retrograde motion. Orbits which visit both domains are intermittently direct or retrograde.
Figures 15(a–b) show two typical (z, pz ) sections of 3D orbits, starting with initial conditions close to two different stable periodic points on the (x, px ) phase planes of the 2D system. In order to obtain the results shown in Figure 15(a), we have taken the point (x 0, p x0) = (0.25, 0) representing approximately the position of the retrograde periodic point in the (x, px ) phase plane of Figure 7(a). Similarly, in Figure 15(b) we can observe the (z, pz ) projections near the retrograde periodic point in the (x, px ) phase plane of Figure 7(b). The position of the periodic point is (x 0, p x0) = (0.25, 0). Note that in both cases the numerical results indicate that for small values of z 0 the motion is regular, while for larger values of z 0 the motion becomes chaotic. The integration time for each ‘invariant curve’ shown in Figures 15(a–b) is 3 × 103 time units. Numerical calculations not given here suggest that the above method can be applied in all regular regions around each retrograde or direct stable periodic point of Figures 7(a–b), 8(a–b), 9(a–b), and 10(a–b).
(a–b) Projections of the sections of 3D orbits with the plane y = 0 when py < 0. The set of the four-dimensional points (x, px , z, pz ) is projected on the (z, pz ) plane.
We must emphasise that the results presented in Figures 15(a–b) are rather qualitative and can be considered as an indication that the transition from regularity to chaos in 3D orbits occurs as the value of z 0 increases. In order to form a more complete and accurate view of the phase space in the 3D system, we computed a large number of 3D orbits (approximately 103) near each periodic point of the (x, px ) phase planes of the 2D system, for different initial conditions (x 0, p x0) and also for different values of z 0. Our target was to determine the average minimum value of z 0 for which the nature of a 3D orbit changes from regular to chaotic. Table 1 shows the value 〈z min〉 near the direct and retrograde stable periodic points of Figures 7(a–b), 8(a–b), 9(a–b), and 10(a–b), for four different values of distance R between the centers of the two galaxies. From Table 1 we can induce three important results. (i) In every case (active or quiet primary galaxy) the minimum value of z 0 in the regions near the retrograde stable periodic points is always larger than the corresponding of the regions of the direct periodic points. (ii) If we compare each kind of galaxy regarding the nuclear activity, (active or quiet) for the same distance R, we observe that when the nucleus is not present at the galactic core (quiet – non-active primary galaxy), the 3D orbits can approach higher values of z 0 and remain regular. On the other hand, when the nucleus is present, the value of 〈z min〉 is smaller. (iii) As the distance between the centers of the galaxies increases, the minimum initial value of z 0 for which a 3D can remain regular increases in both cases (active or quiet primary galaxy). This means that as the two galaxies are at large distances, their mutual interactions are weak enough and, therefore, the majority of the 3D orbits are ordered. Moreover, at large distances, one can conclude that the main factor responsible for the observed chaotic motion is the nuclear activity of the nucleus in the core of the primary galaxy.
Average Value of the Minimum z 0 Near the Direct and Retrograde Periodic Points, for Four Different Values of the Distance R
So far, we have seen that 3D orbits with initial conditions (x 0, p x0, z 0), such that (x 0, p x0) is a point in the chaotic regions of the 2D system, for all permissible values of z 0 are chaotic. On the other hand, the nature (ordered or chaotic) of 3D orbits with initial conditions (x 0, p x0, z 0), such that (x 0, p x0) is a point in the regular regions around the stable direct or retrograde periodic points of the 2D system, depends on the particular value of z 0, as we can see in Table 1. We did not feel that it was necessary to try to define the values of 〈z min〉 for each regular region of the 2D system corresponding to secondary resonances which are represented by sets of multiple small islands of invariant curves in the (x, px ) phase planes shown in Figures 7(a–b), 8(a–b), 9(a–b), and 10(a–b). Numerical results indicate that 3D orbits with initial conditions (x 0, p x0, z 0), such that (x 0, p x0) is a point in the regular regions corresponding to secondary resonances of the 2D system, remain regular for 〈z min〉≲0.052, while for larger values of z 0, they change their character from regular to chaotic.
In Figures 16(a–h) we present eight typical orbits of the 3D dynamical system. Figure 16(a) shows a regular quasi-periodic orbit circulating around the primary galaxy, with initial conditions x 0 = 0.167, y 0 = 0, z 0 = 0.01, and p x0 = p z0 = 0, while the value of p y0 is always found from the Jacobi integral (10). In Figure 16(b), we observe a 3D periodic orbit characteristic of the 1:2 resonance. This orbit has initial conditions x 0 = −0.1895, y 0 = 0, z 0 = 0.02, and p x0 = p z0 = 0. Figure 16(c) shows a 3D periodic orbit with initial conditions x 0 = 0.5788, y 0 = 0, z 0 = 0.01, and p x0 = p z0 = 0, circulating around the center of the primary galaxy. In Figure 16(d), we see a quasi-periodic orbit, with initial conditions x 0 = 0.4295, y 0 = 0, z 0 = 0.015, and p x0 = p z0 = 0. A 3D resonant orbit with initial conditions x 0 = 0.7445, y 0 = 0, z 0 = 0.01, and p x0 = p z0 = 0, moving around the galactic center, is shown in Figure 16(e). In Figure 16(f) we observe a 3D periodic orbit circulating around the primary galaxy. This orbit has initial conditions x 0=−0.9854, y 0=0, z 0=0.01, and p x0=p z0=0. In Figure 16(g) we observe a complicated 3D resonant periodic orbit with initial conditions x 0=0.673, y 0=0, z 0=0.005, and p x0=p z0=0. Figure 16(h) shows a 3D chaotic orbit with initial conditions x 0 = −0.31, y 0 = 0, z 0 = 0.1, and p x0 = p z0 = 0. This orbit goes arbitrarily close to the primary galaxy and it is deflected to higher values of z, on approaching the dense and massive nucleus. We must note that in all 3D orbits shown in Figures 16(a–h), the initial conditions (x 0, p x0) and the values of all the other parameters (M d, M n, R) are as in the corresponding 2D orbits shown in Figures 13(a–h). Moreover, we observe that all regular 3D orbits shown in Figures 16(a– g) stay near the galactic plane and, therefore, support the disk structure of the primary galaxy. The numerical integration time for all 3D orbits shown in Figures 16(a–h) is 200 time units.
(a–h) Eight representative regular orbits of the 3D dynamical system. The initial conditions and more details regarding the values of all the other parameters are given in the text.
5 A THEORETICAL APPROACH
In this section, we shall present some theoretical arguments, together with elementary numerical calculations, in order to explain the numerically obtained relationships given in Figures 11, 12, and 14. Moreover, we will quote different kinds of theoretical techniques in an attempt to explain and study deeper the dynamical structure of the 2D or 3D Hamiltonian system.
Since the potential of the nucleus is integrable with spherical symmetry, the observed chaotic phenomena in our dynamical system should derive mainly from the total |Ft | force. Figure 17 shows a plot of the |Ft | force as a function of the distance R between the two galaxies at the point Pc (x 0, y 0, z 0) = (−0.12, 0.02, 0.01). This point is very close to the center of the primary galaxy, where the chaotic phenomena are more prominent. One may observe from the plot shown in Figure 20 that the two curves, which correspond to two different cases (active and quiet primary galaxy), have almost the same pattern. In both cases, the |Ft | force decreases, tending asymptotically to zero, as the value of the distance R increases. Here, we must point out that the values of the |Ft | force are larger when the primary galaxy is active and, therefore, the corresponding curve is above the one corresponding to the quiet galaxy. Thus, we conclude that when a massive and dense nucleus is present in the center of the primary galaxy, the extent of the chaotic orbits is larger than in the case in which the primary galaxy is quiet. One can see that the pattern between the |Ft | force and the distance R shown in Figure 17 and obtained theoretically is very similar to those connecting the chaotic percentage A% or the LCE with the distance R shown in Figures 11, 12, and 14, respectively, which have been obtained numerically. Thus, the numerically obtained results given in Figures 11, 12, and 14 regarding the evolution of the chaos with respect to the distance R can be justified and explained, in a way, theoretically through the relation between the strength of the total force and the distance R depicted in Figure 17.
A plot of the total |Ft | force as a function of the distance R, when the primary galaxy is active or quiet.
Furthermore, it would be of particular interest to study the structure of the velocity profile that is the plot of the total velocity of a test particle (star),
$\upsilon = \sqrt{\dot{x}^2 + \dot{y}^2 + \dot{z}^2}$
as a function of time, for ordered and chaotic motion. Figure 18(a) shows the velocity profile for a time period of 500 time units corresponding to the regular 3D orbit shown in Figure 16(a). Here, the velocity profile is quasi-periodic and the maximum value of the velocity is about 450 km s−1. This suggests that the regular 3D motion occurs in small velocities. Figure 18(b) shows the velocity profile for a time period of 500 time units corresponding to the chaotic 3D orbit shown in Figure 16(h). In this case, we can discuss two aspects. The first is that the velocity obtains high values up to 790 km s−1 and the second is that the velocity profile appears to be highly asymmetric, displaying abrupt changes and large deviations between the maxima and also between the minima in the [υ−t] plot. Therefore, we may conclude that the chaotic 3D motion occurs in high and abruptly changing velocities.
The physical parameter playing an important role in the orbital behaviour of the stars is the Lz component of the total angular momentum. From our previous experience, we know that on approaching a dense and massive nucleus, low angular momentum stars are scattered off the galactic plane, displaying chaotic motion (Caranicolas & Innanen Reference Caranicolas and Innanen1991; Caranicolas & Papadopoulos Reference Caranicolas and Papadopoulos2003; Caranicolas & Zotos Reference Caranicolas and Zotos2010; Zotos Reference Zotos2011a, Reference Zotos2011c). Of course, here in 3D phase space, things are more complicated than in axially symmetric dynamical models, where the Lz component is conserved. As the motion takes place in a rotating non-axially symmetric system, the Lz component is not conserved and it is given by
\begin{equation}
L_z = x \dot{y} - \dot{x} y - \Omega _p (x^2 + y^2 ).
\end{equation}
\begin{equation}
\langle L_z\rangle = \frac{1}{n} \displaystyle \sum \limits _{i=0}^n L_{zi}.
\end{equation}
In what follows we shall present a semi-theoretical analysis, in order to give a more detailed picture of the structure of the dynamical system and its behaviour. The forces acting on a test particle along the x and y axes are given by the equations
\begin{eqnarray}
F_x &=& -\frac{M_{\rm n} (x - x_1)}{\left(r_{a1}^2 + c_{\rm n}^2\right)^{3/2}} - \frac{M_{\rm d} (x - x_1)}{\left[b^2 + r_{a1}^2 + \left(a + h\right)^2\right]^{3/2}} \nonumber \\
&-& \frac{M_{\rm s} (x - x_2 )}{r_{a2}^3} + \Omega _p^2 x\ - 2\Omega _p \dot{y}, \nonumber\\[-7pt]
\\[-7pt]
F_y &=& -\frac{M_{\rm n} y}{\left(r_{a1}^2 + c_{\rm n}^2\right)^{3/2}} - \frac{M_{\rm d} y}{\left[b^2 + r_{a1}^2 + \left(a + h\right)^2\right]^{3/2}} \nonumber \\
&-& \frac{M_{\rm s} y}{r_{a2}^3} + \Omega _p^2 y + 2\Omega _p \dot{x}. \nonumber
\end{eqnarray}
Contours of the Fx = const (elliptic shaped) together with the contours Fy = const (figure-eight-shaped). Details are given in the text.
At this point, we can expand the semi-theoretical analysis described above, in order to give a more detailed and complete picture of the structure of the 3D dynamical system and its behaviour. The force acting on a test particle along the z axis is given by the equation
\begin{equation}
F_z = - \frac{M_{\rm n} z}{\left(r_1^2 + c_{\rm n}^2\right)^{3/2}} - \frac{M_{\rm d} (a + r_z)z}{r_z \big[b^2 + r_{a1}^2 + (a + r_z )^2\big]^{3/2}} - \frac{M_{\rm s} z}{r_2^3},
\end{equation}
$r_z = \sqrt{h^2 + z^2}$
. It is obvious from Equation (18) that the strength of the Fz
force increases as the mass of the nucleus or the disk increases or their scale lengths decrease. Figure 21(a) shows a 3D plot of the value of the Fz
force on the (x, z) plane. Figure 21(b) depicts the contours of the projections Fz
= const on the (x, z) plane. We observe that for positive values of z, the Fz
force is negative, while for negative values of z, the value of the Fz
is positive, near the center of the primary galaxy. Therefore, lighter colors on the lower half part of the (x, z) plane indicate higher values of the Fz
force, while darker colors on the upper half part of the same plane indicate lower values of the Fz
force. One can observe that near the active galactic core which hosts a dense and massive nucleus, the test particle experiences a very strong Fz
force. The values of all the parameters in Figures 21(a–b) are as in Figure 2.
(a) A 3D plot of the value of the Fz force on the (x, z) plane and (b) contours of the projections Fz = const on the (x, z) plane. Lighter colors indicate higher values of Fz . For positive values of z the Fz force is negative, while for negative values of z the Fz force is positive.
Last but not least, we also investigate the behaviour of the total velocity near the galactic core as a function of the distance x. To do that, we consider the limiting curve, that is, the curve containing all the invariant curves on the (x, px ) phase plane. This can be obtained if we set y = z = py = pz = 0 in Equation (10), yielding
\begin{equation}
\frac{1}{2}p_x^2 = \frac{1}{2}\upsilon ^2 = E_J - \Phi _{\rm t} (x),
\end{equation}
A plot of the total velocity υ(x) as a function of the distance x, for the two cases (active and quiet primary galaxy).
One of the factors responsible for the chaotic motion and other resonance phenomena, such as islandic motion corresponding to sets of multiple islands of invariant curves or invariant tori, are the several inner Lindblad resonances
\begin{equation}
\Omega _p = \Omega - \frac{n}{m} \kappa ,
\end{equation}
Radii of Lindblad Resonances When R = 2.35 and Ω p = 0.419146, for Both Active and Quiet Primary galaxies
Figure 23(a) shows a plot of the curves Ω − nκ/m in the case where the primary galaxy is active, as a function of the radius r, when R = 2.35. We choose this particular value of the distance between the two galaxies because at this distance we can observe a large variety of Lindblad resonances. The numbers 1, 2, 3, and 4 indicate the curves Ω − 2κ/3, Ω − 3κ/4, Ω − 5κ/8, and Ω − 7κ/9, respectively. The straight lines are the curves Ω p = ±0.419146. The values of all the other parameters are as in Figure 2. Figure 23(b) is similar to Figure 23(a) but when the primary galaxy is quiet. Here, the numbers 1, 2, 3, and 4 indicate the curves Ω − 2κ/5, Ω − 3κ/4, Ω − 3κ/7, and Ω − 4κ/9, respectively. As we can see, there are a considerable number of resonance radii for both the direct and retrograde orbits in both cases. Details regarding the resonances and the resonance radii can be obtained from Table 2. In other words, all the Lindblad resonances given in Table 2 are also responsible for the chaotic motion in the primary galaxy. It is also interesting to note that the above resonances produce large chaotic regions for small values of the distance R, while for larger values of R (see Figures 9 a–b) the chaotic regions are extremely small, although the resonance radii are still present. This means that in this case the distance between the two galaxies precedes the Lindblad resonances.
The curves Ω − nκ/m versus r, when R = 2.35 and Ω p = 0.419146, when (a) the primary galaxy is active and (b) the primary galaxy is quiet.
In all cases, the primary galaxy is assumed to have a fixed value of total mass
\begin{equation}
M_{\rm P} = M_{\rm d} + M_{\rm n} = 2.08,
\end{equation}
\begin{eqnarray}
F_x &=& -c_1 -k_1 M_{\rm n}, \nonumber \\
F_y &=& -c_2 -k_2 M_{\rm n}, \nonumber\\[-7pt]
\\[-7pt]
F_z &=& -c_3 -k_3 M_{\rm n}, \nonumber
\end{eqnarray}
\begin{equation}
p_y^2 = 2 [E_J - \Phi _{\rm t}(x) ] - p_x^2 \ge 0.
\end{equation}
\begin{equation}
f (x,p_x) = 2 [E_J - \Phi _{\rm t}(x) ] - p_x^2 = 0
\end{equation}
\begin{equation}
\upsilon = \left[2 E_J + const + c_4 (x_0 - x_1 )^2 + \frac{2 M_{\rm n} (\lambda _1 - \lambda _2 )}{\lambda _1 \lambda _2} \right]^{1/2},
\end{equation}
\begin{eqnarray}
\lambda _1 &=& \sqrt{b^2 + (x_0 - x_1 )^2 + (a + h)^2}, \nonumber\\[-7pt]
\\[-7pt]
\lambda _2 &=& \sqrt{(x_0 - x_1 )^2 + c_{\rm n}^2}. \nonumber
\end{eqnarray}
\begin{equation}
\Theta (M_{\rm n}) = r_0 \sqrt{\frac{M_{\rm n}}{\left(r_0^2 + c_{\rm n}^2\right)^{3/2}} + \frac{M_{\rm d}}{\left(b^2 + r_0^2 + \left(a + h \right)^2\right)^{3/2}}}.
\end{equation}
(a–d) Evolution of different dynamical quantities of the system as a function of the mass of the nucleus M n. Details are given in the text.
6 EVOLUTION OF THE ORBITS IN THE 3D TIME-DEPENDENT MODEL
Let us now follow the evolution of 3D orbits, as mass is transported from the disk of the quiet galaxy to its center. By this procedure, a massive and dense nucleus is developed in the central regions of the primary galaxy. The mass transport is linear following the set of equations:
\begin{eqnarray}
M_{\rm nf}(t) &=& M_{\rm ni} + k t, \nonumber\\[-7pt]
\\[-7pt]
M_{\rm df}(t) &=& M_{\rm di} - k t, \nonumber
\end{eqnarray}
(a) Evolution of the LCE of the 3D orbit in the time-dependent model, following relations (28), and (b) the corresponding velocity profile. The orbit starts as a regular and remains regular during the mass transportation.
(a) Evolution of the LCE of the 3D orbit in the time-dependent model, following relations (28), and (b) the corresponding velocity profile. The orbit starts as a chaotic and remains chaotic during the mass transportation.
(a) Evolution of the LCE of the 3D orbit in the time-dependent model, following relations (28), and (b) the corresponding velocity profile. The orbit starts as a regular but after 100 time units; when the galactic evolution stops it becomes chaotic.
(a) Evolution of the LCE of the 3D orbit in the time-dependent model, following relations (28), and (b) the corresponding velocity profile. The orbit starts as a chaotic but after the galactic evolution, it changes its nature to regular.
In Figure 25(a) we can see the evolution of the LCE for a 3D orbit and for a time period of 105 time units, as the galaxy evolves following the set of equations (28). The initial conditions are x 0 = 0.12, y 0 = 0, z 0 = 0.02, p x0 = 0, and p z0 = 0, while the value of p y0 is found from the Hamiltonian (10) in all cases. The values of all the other parameters are as in Figure 7(b). When t = 100 time units, the mass of the developed nucleus in the core of the galaxy is M nf = 0.08 and the evolution stops. The value of the Hamiltonian is now EJf = −2.60012. The profile of the LCE clearly indicates that this orbit starts as a regular and remains regular during the galactic evolution. Figure 25(b) shows the velocity profile for the same 3D orbit for a time period of 200 time units. The pattern is almost quasi-periodic and therefore the orbit remains regular during the galactic evolution.
On the other hand, things are quite different in Figure 26(a). This orbit has initial conditions x 0 = −0.44, y 0 = 0, z 0 = 0.01, p x0 = 0, and p z0 = 0. The values of all the other parameters are as in Figure 7(b). When t = 100 time units, the mass of the developed nucleus in the core of the primary galaxy is M nf = 0.08 and the mass transportation stops. The value of the Hamiltonian is now EJf = −2.60057. In this case, the profile of the LCE clearly indicates that this orbit starts as a chaotic 3D orbit and remains chaotic during the galactic evolution. Figure 26(b) shows the velocity profile for the same 3D orbit, for a time period of 200 time units. The pattern is highly asymmetric, with a large number of abrupt peaks, large deviations between the maxima and also significant deviations between the minima. Therefore, we conclude that this orbit remains chaotic during the mass transportation.
In Figure 27(a) we observe the evolution of the LCE of a 3D orbit with initial conditions x 0 = −0.196, y 0 = 0, z 0 = 0.01, p x0 = 0, and p z0 = 0. The values of all the other parameters are as in Figure 7(b). After a time interval of 100 time units the mass of the nucleus in the core of the galaxy is M nf = 0.08 and the galactic evolution stops. The Hamiltonian settled to the value EJf = −2.60034. The profile of the LCE given in Figure 27(a) shows that the orbit starts as a regular 3D orbit, but after the galactic evolution it becomes chaotic. It is evident that if mass transport were not present, the orbit would have remained regular. The presence of the nucleus in the core of the primary galaxy has changed the character of the 3D orbit from regular to chaotic. The same result can be obtained by looking at the profile of the velocity of the 3D orbit, as shown in Figure 27(b). For the first 100 time units the pattern is quasi-periodic, indicating regular motion, but after that time interval it becomes highly asymmetric, leading to the conclusion that the orbit has become chaotic.
Figure 28(a) depicts the evolution of the LCE of a 3D orbit with initial conditions x 0 = 0.74, y 0 = 0, z 0 = 0.01, p x0 = 0, and p z0 = 0. The values of all the other parameters are as in Figure 7(b). After a time interval of 100 time units, the mass of the nucleus in the core of the galaxy is M nf = 0.08 and the mass transportation stops. The Hamiltonian settled to the value EJf = −2.60021. The profile of the LCE given in Figure 28(a) indicates that this orbit starts as a chaotic 3D orbit, but after the galactic evolution it becomes an ordered one. Therefore, it is evident that if mass transportation were not present, the orbit this time would have remained chaotic. In this case, the presence of the nucleus in the core of the galaxy has changed the nature of the 3D orbit from chaotic to regular. A similar result can be obtained from the profile of the velocity of the 3D orbit, as shown in Figure 28(b). For the first 100 time units the pattern is highly asymmetric, indicating chaotic motion, but for the rest 100 time units, it becomes almost quasi-periodic, leading to the conclusion that the orbit has now become regular.
The above analysis reveals that as the galaxy tends to create the central nucleus, the character of most 3D regular orbits, starting near the center of the primary galaxy, changes their nature from regular to chaotic. At the same time, the velocity shows significant changes and increases as the mass of the nucleus increases. In other words, as the quiet galaxy shows this active face, the nature of motion presents major alterations. Using the above procedure, we have tested a large number of 3D orbits (approximately 103) in the time-dependent model, describing the formation of a massive and dense nucleus in the core of the primary galaxy, when mass is transported from the disk. Our numerical results indicate that the character of the 3D orbits can change either from regular to chaotic and vice versa or not change at all, as the mass is transported in order to create the nucleus in the central region of the disk galaxy. In particular, from the sample of the 103 tested orbits in the case of the time-dependent model, we conclude that 61% of the orbits altered their nature from regular to chaotic, 22% remain chaotic, 15% remain regular, and only 2% changed their character from chaotic to regular. Thus, it is evident that the formation of a massive nucleus leads the majority of the orbits of the 3D system to become chaotic. Whether the nature of a 3D orbit will change or not during the galactic evolution described by the set of equations (28) depends strongly on the initial conditions (x 0, p x0, z 0) of each orbit. Moreover, as the change of the value of the Hamiltonian (10) is negligible (ΔEJ = |EJf − EJi | ≃ 10−4), we can say that the phase space is transformed to itself during the galactic evolution and, therefore, the orbits can be considered isoenergetic. In order to make this statement more clear, we present an example. If we suppose that the evolving time-dependent model was describing the 2D system, then as the change of the Hamiltonian is negligible, we could say that in the phase plane of the quiet galaxy, shown in Figure 7(b), chaotic regions would appear in the central area and it would be transformed to the phase plane shown in Figure 7(a), as the nucleus is formed in the core of the primary disk galaxy.
7 DISCUSSION AND CONCLUSIONS
In the present article, we have constructed a 3D gravitational model, in order to study the character of motion in a binary galactic system. In particular, we present a dynamical model which is composed of a primary disk galaxy and a small satellite companion. The motion of a test particle (star) in the gravitational field of this galactic system was studied and various techniques were used to identify regular and chaotic motion. The galaxies were assumed to be coplanar and orbiting each other in the same plane on a circular orbit. In particular, this research is indeed an extension to three dimensions of the work of Caranicolas & Innanen (Reference Caranicolas and Innanen2009), who studied a similar two-dimensional binary system of interacting galaxies.
A binary system of interacting galaxies is very complex and, therefore, we need to assume some necessary simplifications and assumptions, in order to be able to study the orbital behaviour of such a complicated stellar system. Thus, 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. Here, we must point out the main restrictions and limitations of our gravitational model: (1) The two galaxies (the primary and the satellite) are assumed to be coplanar, orbiting each other in the same plane on circular orbits. (2) Our dynamical model only deals with the non-dissipative components of the galaxies, stars, or possibly dark matter particles. (3) The potentials we use are rigid and do not respond to the evolving density distribution in a more realistic way. This is because our gravitational model that describes the binary system of the interacting galaxies is not self-consistent. Thus, the two galactic centers remain stationary in the rotating frame chosen. Self-consistent models are usually used when conducting N-body simulations. Obviously, this is out of the scope of the present paper. Once again, note that the above restrictions and limitations of our model are necessary; otherwise it would be extremely difficult, or even impossible, to apply the extensive and detailed dynamical study presented in this study. A similar gravitational model with the same limitations and assumptions was used in Zotos (Reference Zotos2012b) in order to study the motion in a binary system of interacting galaxies, where the second galaxy is no longer treated as a point-mass satellite companion. Thus, we may conclude that the same setup has been successfully applied to describe a much more complex binary galactic system.
In our study, two different cases were investigated: the time-independent model and the time-evolving model, that is the case when mass is transported adiabatically from the disk of the galaxy to its center, forming a massive and dense nucleus. Our numerical calculation indicates that there are several factors responsible for the observed chaotic motion in the time-independent model: (i) the galactic interaction, (ii) the galactic activity, that is the presence of the nucleus, and (iii) the Lindblad resonances. Furthermore, the presence of the nucleus increases the velocities near the central region of the primary galaxy. The value of the velocity depends on the mass of the nucleus and also on the value of its scale length. Regular motion corresponds to low central velocities, while chaotic motion is characterised by high velocities. All the above observations strongly indicate that in the centers of active galaxies, chaotic motion in high velocities is expected. On the other hand, it was found that the two interacting galaxies, for large values of the distance R between their centers, do not present chaotic motion, when the nucleus is not present in the core of the primary galaxy. The results of this work regarding the nature of the orbits and also the factors that affect it or change it are very similar to the corresponding outcomes obtained in Zotos (Reference Zotos2012b) and, therefore, verify in a way our previous conclusions.
We have started our investigation from the Hamiltonian system of two degrees of freedom (2D). Our numerical calculations indicate that in this case, a large part of the phase plane is covered by chaotic orbits, while the regular regions are confined mainly near the stable retrograde periodic point. The chaotic area is larger when the primary galaxy possesses a dense and massive nucleus or when the distance R between the two galaxies is small. It was also found that the velocity near the central region of the galaxy increases significantly when the nucleus is present. For lower values of the total energy and small values of the distance R, no chaos is observed when the nucleus is absent, while a small chaotic region appears in the presence of the nucleus. This means that low-energy stars are in chaotic orbits near the centers of active galaxies with a suitable satellite companion in a circular orbit, while in quiet galaxies they are not. Furthermore, the velocity near the center rises to high values, even for smaller values of the energy. The above discussion strongly indicates that in the centers of active galaxies, chaotic motion and high velocities are expected. This fact, combined with outcomes from previous works, shows that the majority of orbits in galaxies with dense and massive nuclei are in chaotic orbits (see Caranicolas & Innanen Reference Caranicolas and Innanen1991; Caranicolas & Papadopoulos Reference Caranicolas and Papadopoulos2003; Zotos Reference Zotos2012a). This seems reasonable because theoretical results show that the nuclear force near the center increases linearly as M n increases. Furthermore, the velocity near the galactic center strongly depends on the scale size of the nucleus c n. The smaller the c n, the higher is the velocity near the central region of an active galaxy.
To explore and understand the nature of orbits in the 3D dynamical system, we have used our knowledge obtained from the study of the 2D system. Of particular interest was the determination of the regions of initial conditions in the (x, px , z), py < 0, (y = pz = 0) phase space that produce regular or chaotic 3D orbits. As the value of p y0 was found from the Jacobi integral (10), we have used the same value of EJ as in the 2D system and took initial conditions (x 0, p x0, z 0) such that (x 0, p x0) lies in the chaotic regions of the 2D system. It was found that the motion is chaotic for all permissible values of z 0. On the other hand, when (x 0, p x0) was inside a regular region around the direct and retrograde periodic points, the corresponding 3D orbits are regular for small values of z 0, while for larger values of z 0 the orbits become chaotic. The particular values of z 0 for which the transition from regularity to chaos in 3D orbits occurred were different for each regular region of the 2D system. Of particular interest are the results given in Table 1, where we define the average minimum value of z 0 near the direct and retrograde periodic points of the primary galaxy, where the transition from regularity to chaos occurs.
In this research, we have used different kinds of theoretical techniques in an attempt to explain and justify the dynamical structure of the 2D or 3D Hamiltonian system. An important role is played by the Lz component of the test particle's angular momentum. It was found that the values of 〈Lz 〉 for regular 3D orbits are larger than those for chaotic 3D orbits. Thus, the Lz component of the angular momentum is a significant dynamical parameter connected with the regular or chaotic character of orbits in both 2D and 3D dynamical systems. In order to estimate the degree of chaos in the 2D as well as in the 3D dynamical system, we have computed the average value of the LCE for a large number of orbits with different initial conditions in the chaotic regions in each case, for a time period of 105 time units. Our numerical results indicate that the degree of chaos in the 3D binary system of interacting galaxies is smaller than that in the 2D system.
It was of great interest to follow the evolution of orbits as nucleus formed in the central region of the primary disk galaxies. In this procedure, mass was transported from the disk to the nucleus and therefore the quiet galaxy became gradually active, following Equation (28). We observed that the final character of the 3D orbits strongly depends on the particular initial conditions (x 0, p x0, z 0). Thus, regular orbits can change to chaotic and vice versa, or maintain their character (ordered or chaotic) during the nuclear formation and since the nucleus has been formed. It was observed that a number of regular quasi-periodic orbits, starting near the central region of the primary galaxy, become chaotic. This can be seen in a total velocity versus time plot, where the asymmetric profile of the total velocity is in agreement with observational data (see Grosbøl Reference Grosbøl2002), where an increase of the stellar velocity is expected, in regions with significant chaoticity. Moreover, observations show that an asymmetric velocity profile indicates chaotic motion. An interesting question is whether interactions are essential in order to trigger the mass transportation and therefore the galactic activity. This question remains unsolved as astronomers have found binary systems with disk galaxies (see Letawe et al. Reference Letawe2006) showing no signs of interaction at all but harboring active nuclei.
Our theoretical dynamical analysis of the binary galactic system can now be compared with observational data derived from the M51-type binary galactic systems. These systems consist of a large, primary disk galaxy and a smaller satellite companion galaxy moving in a circular orbit around their common mass center (see Klimanov & Reshetnikov Reference Klimanov and Reshetnikov2001). We apply our results to the binary system composed of the primary galaxy NGC 5829 and its small satellite IC 4526. Figure 29 depicts a real image of this binary system. The primary galaxy has an absolute magnitude of MV = −20.66. In order to estimate the disk's radius, we shall use the empirical relation
\begin{equation}
M_V = -6 \log A + 7.14,
\end{equation}
A real image of the binary galactic system composed of the primary galaxy NGC 5829 and its small satellite companion IC 4526.
When studying galactic dynamics, it is invariably found that the stellar rotation velocity remains constant, or flat, with increasing distance away from the galactic center (see Cantinella, Giovanelli, & Haynes Reference Cantinella, Giovanelli and Haynes2006). This result is highly counterintuitive since, based on Newton's law of gravity, the rotational velocity would steadily decrease for stars farther away from the galactic center. By this particular argument, the flat rotational curves seem to imply that each galaxy must be surrounded by significant amounts of dark matter. It has been postulated and generally accepted that the dark matter would have to be located in a massive halo enshrouding each galaxy (see Caranicolas & Zotos Reference Caranicolas and Zotos2009b, Reference Caranicolas and Zotos2011a, Reference Caranicolas and Zotos2011b). The first real surprise in the study of dark matter lay in the outermost parts of galaxies known as galaxy halos. Here, there is negligible luminosity, yet there are occasional orbiting gas clouds, which allow us to measure rotation velocities and distances. The rotation velocity was found not to decrease with increasing distance from the galactic center, implying that the mass distribution of the galaxy cannot be considered like the light distribution. The mass must continue to increase since the rotation velocity satisfies
\begin{equation}
\Theta ^2 = \frac{G M}{r},
\end{equation}
Forty years ago, galactic activity and interactions between galaxies were viewed as unusual and rare. Nowadays, they seem to be segments in the life of many galaxies. From the astrophysical point of view, in the present work, we have tried to connect galactic activity and galactic interactions with the nature of orbits (regular or chaotic) and also with the behaviour of the velocities of stars in the primary galaxy. We consider the outcomes of the present research to be an initial effort, in order to explore the dynamical structure of the 3D binary stellar system in more detail. As results are positive, further investigation will be initiated to study all the available phase spaces, including orbital eccentricity (elliptic orbits) of the small companion galaxy and its inclinations to the primary galaxy. Moreover, we shall try to use the outcomes obtained from this initial and simple dynamical model to conduct computer N-body simulations in a binary system of interacting galaxies (a primary galaxy with a satellite companion), in order to reveal changes in their orbital properties through merger processes and tidal effects, which are obviously out of the scope of this study.
ACKNOWLEDGMENTS
The author would like to express his warmest thanks to the anonymous referee for careful reading of the manuscript and for his very useful and illuminating suggestions and comments that greatly improved the quality and also the clarity of the present article.
