Orbital dynamics in realistic galaxy models: NGC 3726, NGC 3877 and NGC 4010

TL;DR

This study models the orbital dynamics of stars in three spiral galaxies using a generalized Miyamoto-Nagai potential, revealing conditions for regular and chaotic stellar orbits based on angular momentum and energy.

Contribution

It introduces a generalized potential model tailored to specific galaxies and analyzes their orbital dynamics, highlighting the presence of chaotic stellar orbits.

Findings

Regular orbits for high angular momentum or low energy

Chaotic orbits for low angular momentum or high energy

Model enables identification of chaotic bounded stellar orbits

Abstract

In the present paper, using a generalization of the Miyamoto and Nagai potential we adjusted the observed rotation curves of three specific spiral galaxies to the analytical circular velocities. The observational data have been taken from a 21 cm-line synthesis imaging survey using the Westerbork Synthesis Radio Telescope, for three particular galaxies in the Ursa Major cluster: NGC 3726, NGC 3877 and NGC 4010. Accordingly, the dynamics of the system is analyzed in terms of the Poincar\'e sections method, finding that for larger values of the angular momentum of the test particle or lower values its total energy the dynamics is mainly regular, while on the opposite cases, the dynamics is mainly chaotic. Our toy model opens the possibility to find chaotic bounded orbits for stars in those particular galaxies.

Equations

$${\nabla }^2\Phi (r,\theta )=\frac{1}{r^2}\frac{\partial }{\partial r}\left(r^2\frac{\partial \Phi }{\partial r}\right)+\frac{1}{r^2\sin \theta }\frac{\partial }{\partial \theta }\left(\sin \theta \frac{\partial \Phi }{\partial \theta }\right)=0,$$

$$\Phi (r,\theta )=\sum _{l=0}^{\infty }\left(A_l{r}^l-B_l{r}^{-(l+1)}\right)P_l(\cos \theta ),$$

$$\Phi (r,\theta )=-\sum _{l=0}^{\infty }\frac{B_l{P}_l(\cos \theta )}{r^{l+1}}.$$

$$\cos \theta =z/r\mathrm{and}r=\sqrt{R^2+z^2},$$

$$z\to z^{\ast }=a+\sqrt{z^2+b^2},$$

$$\Phi (R,z)$$

$$\Sigma =\frac{1}{4\pi G}\left(\frac{{\partial }^2\Phi }{\partial R^2}+\frac{1}{R}\frac{\partial \Phi }{\partial R}+\frac{{\partial }^2\Phi }{\partial z^2}\right),$$

$$v^2$$

$$k^2$$

$${\nu }^2$$

$$\mathcal{L}=\frac{1}{2}\left[{\dot{R}}^2+(R\dot{\varphi }{)}^2+{\dot{z}}^2\right]-\Phi (R,z),$$

$$p_R=\dot{R},p_{\varphi }=R^2\dot{\varphi },p_z=\dot{z},$$

$$\mathcal{H}=\frac{1}{2}\left(p_R^2+p_z^2\right)+{\Phi }_{\mathrm{eff}}(R,z),$$

$${\Phi }_{\mathrm{eff}}(R,z)=\frac{L_z^2}{2R^2}+\Phi (R,z).$$

$$\dot{R}$$

$$\dot{z}$$

$$\dot{p_R}$$

$$\dot{p_z}$$

$$\mathcal{H}(R,z,p_R,p_z)=\mathcal{H}(R_0,z_0,p_{R_0},p_{z_0})=h,$$

Full text

Physical Sciences

**Orbital dynamics in realistic galaxy models: NGC 3726, NGC 3877 and NGC 4010

**

Fredy L. Dubeibe1, Sandra M. Martínez-Sicachá2, Guillermo A. González2111Correspondencia: G. A. González, guillermo.gonzalez@saber.uis.edu.co, Recibido: 5 de octubre de 2018; Aceptado: 14 de febrero de 2019.

1Grupo de Investigación Cavendish, Facultad de Ciencias Humanas y de la Educación, Universidad de los Llanos, Villavicencio, Colombia

2Grupo de Investigación en Relatividad y Gravitación, Escuela de Física, Universidad Industrial de Santander, A.A. 678, Bucaramanga, Colombia

Abstract

In the present paper, using a generalization of the Miyamoto and Nagai potential we adjusted the observed rotation curves of three specific spiral galaxies to the analytical circular velocities. The observational data have been taken from a 21 cm-line synthesis imaging survey using the Westerbork Synthesis Radio Telescope, for three particular galaxies in the Ursa Major cluster: NGC 3726, NGC 3877 and NGC 4010. Accordingly, the dynamics of the system is analyzed in terms of the Poincaré sections method, finding that for larger values of the angular momentum of the test particle or lower values its total energy the dynamics is mainly regular, while on the opposite cases, the dynamics is mainly chaotic. Our toy model opens the possibility to find chaotic bounded orbits for stars in those particular galaxies.

Key words: Stellar dynamics; Galaxies: kinematics and dynamics; Nonlinear dynamics and chaos.

Dinámica orbital en modelos realistas de galáxias: NGC 3726, NGC 3877 y NGC 4010

Resumen

En el presente trabajo, utilizando una generalización del potencial de Miyamoto-Nagai, se ajustan las curvas de rotación observadas de tres galaxias espirales a las velocidades circulares analíticas. Los datos observacionales se tomaron de un conjunto de imágenes de línea de 21 centímetros (o línea HI) obtenidos con el Westerbork Synthesis Radio Telescope (WSRT), para tres galaxias particulares en el grupo de la Ursa Major: NGC 3726, NGC 3877 y NGC 4010. Seguidamente, se analiza la dinámica del sistema en términos del método de secciones de Poincaré, encontrando que para valores grandes del momento angular de la partícula de prueba o valores bajos su energía total, la dinámica es principalmente regular, mientras que en los casos opuestos, la dinámica es principalmente caótica. Nuestro modelo abre la posibilidad de encontrar órbitas caóticas acotadas para estrellas presentes en esas galaxias particulares.

Palabras clave: Dinámica estelar, Galaxias: cinemática y dinámica, Dinámica no lineal y caos.

Introduction

Since the seminal paper by Miyamoto & Nagai (1975), the literature on three-dimensional analytical models for the gravitational field of different types of galaxies has grown considerably. In this respect, particular attention deserve the models proposed by Jaffe (1983) and Hernquist (1990), who derived analytical models which closely approximate the light distribution for spherical and elliptical galaxies, respectively. A few years later, Long & Murali (1992) presented an analytical potential for barred galaxies that reduces to the Miyamoto-Nagai disk by an appropriate setting of the free parameters, while Dehnen (1993) generalized the Jaffe and Hernquist models by means of a family of density profiles with different central slopes. More recently, Vogt & Letelier (2005) derived an analytical expression for the gravitational field of galaxies, based on the multipolar expansion up to the quadrupole term. Using a different approach, **González **et al. (2010) obtained a family of finite thin-discs models for four galaxies in the Ursa major cluster in which the circular velocities were adjusted to fit the observed rotation curves.

One advantage of an analytical galaxy model is the possibility to study the dynamics (regular or chaotic) of orbits. This can be considered one of the standing problems in galactic dynamics because it could allow us to understand the formation and evolution of galaxies (Contopoulos, 1979), as shown by the pioneer simulations of Lindblad (1960). Despite the fact that early papers on this topic studied only regular orbits in the meridional plane (Martinet & Mayer, 1975, Manabe, 1979, Greiner, 1987, Lees & Schwarzschild, 1992), soon after, the existence of chaos on the orbital motion started to be considered by Caranicolas (1996) and Caranicolas & Papadopoulos (2003). In the majority of cases all these studies focused on the distinction between regular and chaotic orbits (Manos & Athanassoula, 2011, Bountis ****et. al., 2012, **Manos **et al., 2013) or the influence of the galaxy components (nucleus, bulge, disk, halo) on the character of orbits, see e.g. (Zotos, 2012, Zotos & Caranicolas, 2013, Zotos, 2014). Notwithstanding the evidence that both chaotic and regular motions are possible in many axisymmetric potentials, recent studies on generalized axisymmetric potentials suggest that a third integral of motion seems to exist for energy values closer to the escape energy (**Dubeibe **et al., 2018, **Zotos **et al., 2018). Hence, such apparent ambiguity might only be solved by performing systematic studies of each particular model.

In this paper, we are interested in meridional motions of free test particles (stars) in presence of analytical realistic galaxy models. Our models possess axial symmetry, which is a good approximation given the morphology of galaxies that are mainly approximate figures of revolution. Additionally, the galaxy components were not added one by one, instead of this, we derived a generalized Miyamoto-Nagai model that can be adjusted very accurately to fit the observed rotation curve and hence it is assumed that all (or most of) the components are taken into account. The determination of the specific values of the coefficients of the series expansion let us calculate the corresponding surface densities and all the kinematic quantities characterizing the particular galaxy models. Unlike the models derived by **González **et al. (2010), which exhibit instabilities to small vertical perturbations (see e.g. the cases of NGC 3877 and NGC 4010), our models satisfy the stability conditions for radial and vertical perturbations. On the other hand, the dynamics of the orbits is studied through the Poincaré surfaces of section, showing that the orbital motion exhibits a strong dependence on the angular momentum and energy of the test particles (stars).

The paper is organized as follows: in the first section, we derive the generalized Miyamoto-Nagai model. Next, from the new potential the explicit expressions for the physical quantities of interest are determined. In the second section we adjust the observed rotation curves of three specific spiral galaxies (NGC 3726, NGC 3877 and NGC 4010) to the analytical circular velocities derived with our model. Then, the mass-density profiles are calculated, along with the vertical and epicyclic frequencies, showing that our model not only is well-behaved but also satisfy the stability conditions. A dynamical analysis in terms of the Poincaré surfaces of section is performed in the third section. Finally, in the fourth section, we summarize our main conclusions.

Generalized Miyamoto-Nagai model

Let us start considering the axially symmetric Laplace’s equation in spherical coordinates

[TABLE]

whose general solution reads as

[TABLE]

where $A_{l}$ and $B_{l}$ are constants to be determined, $P_{l}$ are the Legendre polynomials, and the notation $(r,\theta,\phi)$ means (radial, polar, azimuthal) coordinates, respectively.

Since $\Phi(r,\theta)$ denotes the gravitational potential of an axisymmetric finite distribution of mass, the boundary condition $\lim_{r\to\infty}\Phi(r,\theta)=0$ must be satisfied, thus the solution (2) takes the form

[TABLE]

Following Vogt & Letelier (2005), in order to obtain a generalized Miyamoto-Nagai model and for the sake of simplicity, we shall consider terms up to $l=3$ in (3), therefore, transforming to cylindrical coordinates $(R,z)$ by means of the relations

[TABLE]

and applying the additional transformation (Satoh, 1980),

[TABLE]

with $a$ and $b$ two arbitrary parameters, the generalized potential takes the form111It should be noted that setting $B_{0}=GM,B_{1}=B_{2}=B_{3}=0$ in (6), we get the well-known Miyamoto-Nagai Potential (Miyamoto & Nagai, 1975).

[TABLE]

Once the potential has been specified, the mass-density distribution $\Sigma$ can be calculated directly from Poisson equation,

[TABLE]

while the circular velocity $v$ of particles in the galactic plane, the epicyclic frequency $k$, and the vertical frequency $\nu$ of small oscillations about the equilibrium circular orbit, can be obtained from the following expressions evaluated at $z=0$ (Binney & Tremaine, 2011)

[TABLE]

From (8-10), it is important to emphasize that a feasible model must satisfy the constraints set by the conditions $v^{2}\geq 0$, $k^{2}\geq 0$, and $\nu^{2}\geq 0$, where the last two inequalities are understood as stability conditions (Vogt & Letelier, 2005).

As is evident from the preceding paragraphs, the galactic models and its associated physical quantities are uniquely determined by the set of constants $a,b,B_{0},B_{1},B_{2}$, and $B_{3}$, which (taking a pragmatic approach) can be estimated from the observational data of the corresponding rotation curves, as we will discuss in detail in the next section.

Rotation curves fitting

The observational data were taken from Verheijen & Sancisi (2001) for three specific galaxies in the Ursa Major cluster: NGC 3726, NGC 3877 and NGC 4010. Following the procedure outlined in **González **et al. (2010), we take the galaxy radius $R_{d}$ as the given by the largest tabulated value of the data. Thus, introducing dimensionless variables $\tilde{R}=R/R_{d},\tilde{z}=z/R_{d},\tilde{a}=a/R_{d},\tilde{b}=b/R_{d}$ and setting $\tilde{B_{0}}=B_{0}/R_{d},\tilde{B_{1}}=B_{1}/R_{d}^{2},\tilde{B_{2}}=B_{2}/R_{d}^{3}$, and $\tilde{B_{3}}=B_{3}/R_{d}^{4}$, the nonlinear least square curve fitting method allows us to calculate the numerical values of the parameters for each particular galaxy. The resulting values of $\tilde{a},\tilde{b},\tilde{B_{0}},\tilde{B_{1}},\tilde{B_{2}}$, and $\tilde{B_{3}}$, for the three galaxies under consideration, are given in Table 1.

[TABLE]

Table 1. Parameters for each particular galaxy model.

Figure 1. Model fitted to the rotation curve of NGC 3726 using the parameters given in the first column of Table 1. (a) The solid curve indicates the rotation velocity calculated from (8) while the error bars denote the velocity dispersions of the observational data. (b) Normalized mass-density distribution $\tilde{\Sigma}$ at $z=0$, calculated from (7). (c) Constant-density curves of equation (7) in the meridional plane. (d) Epicyclic frequency (9) evaluated on $z=0$. (e) Vertical frequency (10) evaluated on $z=0$.

Figure 2. Model fitted to the rotation curve of NGC 3877 using the parameters given in the second column of Table 1. (a) The solid curve indicates the rotation velocity calculated from (8) while the error bars denote the velocity dispersions of the observational data. (b) Normalized mass-density distribution $\tilde{\Sigma}$ at $z=0$, calculated from (7). (c) Constant-density curves of equation (7) in the meridional plane. (d) Epicyclic frequency (9) evaluated on $z=0$. (e) Vertical frequency (10) evaluated on $z=0$.

In panels (a) of Figures 1, 2, and 3, we show the observational data (points) of the rotation curve with the corresponding velocity dispersions (error bars) as reported by Verheijen & Sancisi (2001) for NGC 3726, NGC 3877 and NGC 4010. The solid lines correspond to the analytical expressions (8) fitted to the rotation curves. As can be seen, in each case the model fits the observed data with good accuracy. Additionally, in panels (b) of Figures 1, 2, and 3, we plot the normalized mass-density distribution (7) at $z=0$ for the three galaxies, as a function of the dimensionless radial coordinate $\tilde{R}$. Here, we obtain a well-behaved mass-density function, showing a maximum value at the center that decreases to zero at the edge of the disk. On the other hand, in panels (c) of Figures 1, 2, and 3, we present four isodensity curves of the mass-density distribution (7) in the meridional plane ($\tilde{R},\tilde{z}$), showing that each model corresponds to a very different mass distribution. Finally, from panels (d) and (e) of the same figures, it is noteworthy that in the three cases the stability conditions are fully satisfied.

Figure 3. Model fitted to the rotation curve of NGC 4010 using the parameters given in the third column of Table 1. (a) The solid curve indicates the rotation velocity calculated from (8) while the error bars denote the velocity dispersions of the observational data. (b) Normalized mass-density distribution $\tilde{\Sigma}$ at $z=0$, calculated from (7). (c) Constant-density curves of equation (7) in the meridional plane. (d) Epicyclic frequency (9) evaluated on $z=0$. (e) Vertical frequency (10) evaluated on $z=0$.

Stellar Dynamics

It is a well-known fact that using rough estimates of the dimensions of typical stars and galaxies, the collision interval between stars is about $10^{8}$ times longer than the average age for most galaxies (Binney & Tremaine, 2011). This implies that the star’s motion can be determined solely by the gravitational attraction of the galaxy and that collisions between stars are so rare that are irrelevant (Maoz, 2016). Therefore, as a first approximation, the orbital dynamics of a star in a given galaxy can be studied following the usual Lagrangian and Hamiltonian approaches for the motion of a test particle in the presence of an estimated gravitational potential.

The orbital motion of a test particle in an axisymmetric potential is governed by the Lagrangian

[TABLE]

with $(R,\phi,z)$ the usual cylindrical coordinates. The generalized canonical momenta read as

[TABLE]

and the Hamiltonian takes the form

[TABLE]

with

[TABLE]

Here, $L_{z}=p_{\phi}=$constant, denotes the conserved component of angular momentum about the $z$-axis.

From (13), the resulting Hamilton’s equations of motion can be expressed as

[TABLE]

where $\Phi(R,z)$ is given by Eq. (6) and its respective parameters should be taken from Table 1.

Since the Hamiltonian is autonomous, $\mathcal{H}$ is an integral of motion

[TABLE]

with $h$ the energy of an orbit.

The existence of an analytic integral of motion reduces the phase space dimensionality, and hence the Poincaré surface of section is an appropriate and well-established method to analyze the dynamics of the system. Taking into account the axial symmetry associated to the system, it is customary to choose the equatorial plane $z=0$ as the Poincaré plane in order to represent the surface of sections in the $(\tilde{R},\dot{\tilde{R}})$-plane. The orbits were numerically integrated forward in time for 1000 units of time by using a Runge-Kutta-Fehlberg Method (RKF45), with this setting the numerical error related to the conservation of the energy is at most $10^{-14}$. In all cases we set $z_{0}=p_{R_{0}}=0$ and we scan the phase space with a large number of initial conditions for the radii $R_{0}$, these three values allow us to determine the values of $p_{z_{0}}$ through the relation (19).

Figure 4. Poincaré surfaces of section of NGC 3726 for different values of angular momentum $L_{z}$ with $h=-1$.

Figure 5. Poincaré surfaces of section of NGC 3726 for different values of total energy $h$ with $L_{z}=1$.

Figure 6. Poincaré surfaces of section of NGC 3877 for different values of angular momentum $L_{z}$ with $h=-1$.

Figure 7. Poincaré surfaces of section of NGC 3877 for different values of total energy $h$ with $L_{z}=1$.

Figure 8. Poincaré surfaces of section of NGC 4010 for different values of angular momentum $L_{z}$ with $h=-10$.

Figure 9. Poincaré surfaces of section of NGC 4010 for different values of total energy $h$ with $L_{z}=5$.

The transition from regularity to chaos (or viceversa) that takes place for the three considered galaxy models was inspected through the Poincaré sections in Figs. 4-9, by using different values of $L_{z}$ (Figs. 4, 6 and 8) and $h$ (Figs. 5, 7 and 9). It can be observed that the orbital motion exhibits a strong dependence on the angular momentum $L_{z}$ and energy $h$ of the test particle. In particular, from the surfaces of section presented in Figs. 4, 6, and 8, we may infer that there exists an increase in the regularity of the system for larger values of the angular momentum $L_{z}$, i.e. if there exists a chaotic sea the increase of $L_{z}$ will fill the phase with KAM islands, while the opposite effect is observed for larger values of energy $h$ (see Figs. 5, 7, and 9), where the KAM islands deform and shrink giving place to larger regions of chaos.

Concluding remarks

In the present paper, using the general solution to the Laplace equation, we have derived a generalized Miyamoto-Nagai potential. By means of the nonlinear least square fitting, the analytical velocity curves were adjusted to the observed ones of three specific spiral galaxies: NGC 3726, NGC 3877 and NGC 4010. The resulting analytical models were used to determine the mass-density distributions and the vertical and epicyclic frequencies, showing that unlike the results presented in **González **et al. (2010) for NGC 3877 and NGC 4010, our models satisfy the stability conditions for radial and vertical perturbations. Even though the set of models presented here should be considered as a rough approximation, the circular velocities were shown to fit very accurately to the observed rotation curves and in the three cases the stability conditions were fully satisfied. Here, it is important to note that contrary to the observed tendency in the Miyamoto-Nagay model, where the limit $a\rightarrow 0$ reduces to the Plummer sphere, our models exhibit a tendency to an spherical mass distribution with increasing of the parameter $a$.

On the other hand, by using the Poincaré section method we have also studied the dynamics of the meridional orbits of stars in presence of the gravitational field of the galaxy models. From our results it may be inferred that there exists an increase in the regularity of the orbits for larger values of the angular momentum, while for larger values of energy the orbits tend to be more chaotic. Our toy models suggest that in the three galaxy models chaotic orbits are possible, however the chaotic behavior is very weak for the NGC 3877 model in comparison to NGC 3726 and NGC 4010. It should be noted that none of the studied models showed a fully chaotic phase space. Our results could have significant implications for the study of the dynamics and kinematics of these three specific galaxies, since the regular or chaotic behaviors could shed lights into the evolution and structure of these galaxies, i.e., in phase space, regular orbits are trapped in the vicinity of neighbor orbits, while chaotic orbits, by its own nature, will diverge exponentially in time from its neighbors by filling the phase space in an erratic manner.

Acknowledgments

We would like to thank the anonymous referees for their useful comments and remarks, which improved the clarity and quality of the manuscript. FLD, SMM and GAG gratefully acknowledges the financial support provided by COLCIENCIAS (Colombia) under Grants No. 8840 and 8863.

Authors’ contributions

All authors make substantial contributions to conception, design, analysis and interpretation of data. All authors participate in drafting the article and reviewed the final manuscript.

Conflict of interest

The authors declare that they have no conflict of interest.