# Shape and orientation of stellar velocity ellipsoids in spiral galaxies

###### Abstract

We present a numerical study of the properties of the stellar velocity distribution in stellar discs which have developed a saturated, two-armed spiral structure. We follow the growth of the spiral structure deeply into the non-linear regime by solving the Boltzmann moment equations up to second order. By adopting the thin-disc approximation, we restrict our study of the stellar velocity distribution to the plane of the stellar disc. We find that the outer (convex) edges of stellar spiral arms are characterized by peculiar properties of the stellar velocity ellipsoids, which make them distinct from most other galactic regions. In particular, the ratio of the smallest versus largest principal axes of the stellar velocity ellipsoid can become abnormally small (as compared to the rest of the disc) near the outer edges of spiral arms. Moreover, the epicycle approximation fails to reproduce the ratio of the tangential versus radial velocity dispersions in these regions. These peculiar properties of the stellar velocity distribution are caused by large-scale non-circular motions of stars, which in turn are triggered by the non-axisymmetric gravitational field of stellar spiral arms.

The magnitude of the vertex deviation appears to correlate globally with the amplitude of the spiral stellar density perturbations. However, locally there is no simple correlation between the vertex deviation and the density perturbations. In particular, the local vertex deviation does not correlate with the local gravitational potential and shows only a weak correlation with the local gravitational potential gradient. We find that the local vertex deviation correlates best with the spatial gradients of mean stellar velocities, in particular with the radial gradient of the mean radial velocity.

###### keywords:

Galaxies: general – galaxies: evolution – galaxies: spiral – galaxies: kinematics and dynamics## 1 Introduction

Already since the end of the nineteenth century studies of stellar kinematics in the solar neighbourhood revealed that the velocity distribution of stars in galactic discs is non-isotropic (Kobold kobold90 ). An early interpretation by Kapteyn (kapteyn05 ) and Eddington (eddington06 ) invoked a superposition of (isotropic) stellar streams with different mean velocities, by this creating an anisotropic velocity distribution. However, an alternative interpretation by Schwarzschild (schwarzschild07 ) became the general framework for describing the local velocity distribution of stars of equal age, though a contamination by stellar streams is still considered to be important for a proper analysis of velocity data. The Schwarzschild distribution is based on a single but anisotropic ellipsoidal distribution function. It is characterized by a Gaussian distribution in all three directions (radial ), (tangential ) and (vertical ) in velocity space, but with different velocity dispersions , , and . In general, the velocity distribution is described by a velocity dispersion tensor

(1) |

where and denote the different coordinate directions, i.e. radial (), azimuthal () and vertical () direction
and gives the corresponding velocity. The bar denotes local
averaging over velocity space^{1}^{1}1We emphasize that
our definition of the velocity dispersion tensor is
different from the definition of e.g. Kuijken & Tremaine (KT )
who do not use a square mark. Hence, our is equivalent to
Kuijken & Tremaine’s . .
The principal axes of this tensor form an imaginary surface that is called the
velocity ellipsoid. The principal axes of the velocity ellipsoid need not to be aligned with
the coordinate axes. In that case, the vertex deviation (which is
defined as the angle between the direction in velocity space
pointing from the Sun to the Galactic centre and
the direction of the major principal axis) does not
vanish.

The general trend for the velocity dispersions is that the radial velocity dispersions are the largest and the vertical velocity dispersions are the smallest among stars of same age and stellar type. This trend holds when dispersions are averaged locally, i.e. in the solar neighbourhood (Dehnen & Binney Dehnen ), or globally over the entire galactic discs. For instance, the ratio of the radial, tangential, and vertical stellar velocity dispersions in the disc of NGC 488 is approximately (Gerssen et al. Gerssen ).

In the case of stationary, axisymmetric systems and appropriate
distribution functions^{2}^{2}2i.e. those DFs only depending on the
two integrals of motion energy and the -component of
the angular momentum or those obeying a third integral of motion
which is symmetric under a simultaneous change of sign of and
(DF) the velocity ellipsoid is aligned with the coordinate
axes. If the stellar orbits are nearly circular (i.e. with small
radial amplitudes with respect to the guiding centre), then the Oort ratio of tangential versus radial
stellar velocity dispersions can be calculated in zero order
approximation within the epicyclic theory to be , where and
are the usual Oort constants. For a flat rotation curve one gets
. More realistic
axisymmetric distribution functions give first order corrections,
which yield values of (or 0.66) for the solar vicinity (Kuijken
& Tremaine kuijken91 , hereafter KT91). Clearly, these values
exceed those predicted by the epicycle approximation. On the other hand,
observations in
the Milky Way (MW) give lower (mean) values around , which are
substantially below the epicycle value (Kerr & Lynden-Bell
kerr86 ; Ratnatunga & Upgren ratnatunga97 ). The exact
reason for this discrepancy is still unknown. It might stem from a DF
asymmetric in (KT91) or from invalid assumptions, i.e. a
violation of stationarity or axisymmetry.

Additionally, non-vanishing vertex deviations were found in the solar vicinity in many studies (e.g. Wielen wielen74 , Mayor mayor72 , Ratnatunga & Upgren ratnatunga97 ). These vertex deviations could be explained by spiral structures (Mayor mayor70 , Yuan yuan71 ), by this supporting the importance of non-axisymmetry of the underlying gravitational potential. Moreover, a clear relation between and has been found. For smaller velocity dispersions larger vertex deviations have been measured, e.g. for a radial velocity dispersion of about 15 km s one gets a vertex deviation of (Wielen wielen74 ).

A major step in analysing the local stellar velocity space was the astrometric satellite mission Hipparcos (ESA esa97 ). Earlier results concerning the general behaviour of the velocity ellipsoid, e.g. the dependence of the dispersions or the vertex deviations as a function of were corroborated (Dehnen & Binney Dehnen , Bienaymé bienayme99 , Hogg et al. hogg05 ). Moreover, the more numerous known proper motions allowed for a detailed mapping of the velocity space in the solar vicinity. It turned out that the local velocity distribution of stars exhibits a lot of substructure (Dehnen dehnen98 , Alcobé & Cubarsi alcobe05 ). For example, two major peaks in the velocity distribution were found where the smaller secondary peak is well detached by at least 30 km s from the main peak (the “u anomaly”). Dehnen (dehnen00 ) attributed this bimodality to the perturbation exerted by the Galactic bar assuming that its outer Lindblad resonance (OLR) is close to the Sun. Mühlbauer & Dehnen (muehlbauer03 ) showed that a central bar can also explain vertex deviations of about 10, and Oort ratios below 1/2. Recent investigations about the influence of the Galactic bar on Oort’s constant corroborated these results and gave stronger limits on the bar’s properties (Minchev et al. minchev07b ). Minchev & Quillen (minchev07a ) showed also that a stationary spiral cannot reproduce the ”u anomaly”. Though less likely, an alternative interpretation of the bimodality as a result of non-stationary spiral arms cannot be ruled out completely.

In general, a non-axisymmetric gravitational potential might lead to both, a misalignment of the velocity ellipsoid and an additional correction of the velocity dispersion ratio with respect to the value predicted by the epicycle approximation. In the case of spiral perturbations, this conclusion was made e.g. by Mayor (mayor70 , for a review of different mechanisms see KT91) and numerically confirmed recently by Vorobyov & Theis (VT , hereafter VT06). Using linear perturbation analysis KT91 studied different discs subject to oscillating perturbations (like a rotating bar or a spiral) or subject to a static non-axisymmetric dark matter halo. In the case of large-scale oscillating perturbations the predicted vertex deviations oscillate, too. Once the axisymmetric background gravitational potential and the pattern speed of the perturbation are fixed, its amplitude depends in a complex way only on the perturbed mean velocity and its spatial gradient. For a static elliptic potential Kuijken & Tremaine (KT ) showed that the observed mean kinematic data ( and ) can be explained, if the ellipticity of the disc is about and the Sun’s location is near to the minor axis of the corresponding potential. It is remarkable that the correction terms introduced by the non-axisymmetry are negative (i.e. in agreement with the observations in the Milky Way), whereas the higher order terms derived from a more realistic axisymmetric DF are positive. Similarly, Blitz & Spergel (blitz91 ) argued for a slowly rotating mildly triaxial Galactic potential, resulting in an outward motion of 14 km s of the LSR and a vertex deviation of 10.

However, it is not clear if the non-axisymmetric gravitational field is entirely responsible for the observed vertex deviation. The existence of moving groups of stars was shown to produce large vertex deviations (Binney & Merrifield BM ). The situation may become even more complicated because moving groups of stars may in turn be caused by the non-axisymmetric gravitational field of spiral arms. Therefore, a detailed numerical study of the vertex deviation in spiral galaxies is necessary.

In this paper, we present the stellar hydrodynamics simulations of linear and non-linear stages of the spiral structure formation in stellar discs by solving numerically the Boltzmann moment equations up to second order in the thin-disc approximation. In contrast to linear stability analyses, we follow the growth of the spiral structure self-consistently. Compared to -body simulations, we have the advantage that we can follow the spiral structure formation from a very low perturbation level deeply into the non-linear regime. The former is usually hardly accessible by self-consistent -body simulations due to particle noise. We study the local evolution of the velocity dispersion tensor in the plane of the disc, which includes the anisotropy (), the ratio of the smallest versus largest principal axes of the stellar velocity ellipsoid, and the vertex deviation ().

The paper is organized as follows. In Section 2 we give a short description of the numerical code. The model galaxy and initial conditions are described in Section 3. The physical explanation for the growth of the spiral structure in our model disc is given in Section 4. The shape and orientation of stellar velocity ellipsoids is studied in Section 5. In Section 6 we discuss the applicability of the epicycle approximation and possible causes of the vertex deviation in spiral galaxies. The summary and conclusions are given in Section 7.

## 2 The method and numerical code

Stars are collisionless objects that move on orbits determined by the large-scale gravitational potential. An exact description of such a system requires the solution of the collisionless Boltzmann equation for the distribution function of stars in the phase space (). A general solution of the Boltzmann equation, which is defined in a six-dimensional position and momentum phase space, is usually prohibited due to the large computational expense. However, taking moments of the Boltzmann equation in velocity space turns out to yield a set of equations that is numerically tractable.

In our recent paper (VT06) we presented the BEADS-2D code that is designed to study the dynamics of stellar discs in galaxies. The BEADS-2D code is a finite-difference numerical code that solves the Boltzmann moment equations up to second order in the thin-disc approximation on a polar grid (). More specifically, the BEADS-2D code solves for the stellar surface density , mean radial and tangential stellar velocities and , and stellar velocity dispersion tensor. The latter includes the squared radial and tangential stellar velocity dispersions and , respectively, and the mixed velocity dispersion . We close the system of Boltzmann moment equations by adopting the zero-heat-flux approximation. This approximation is equivalent to assuming the absence of any heat transfer in the classic fluid dynamics approach. The BEADS-2D code has shown a good agreement with the predictions of linear stability analysis of self-gravitating stellar discs (e.g. Polyachenko et al. Polyachenko , Bertin et al. bertin00 ).

The use of the Boltzmann moment equation approach allows us to study non-isotropic stellar systems, in particular, the shape and orientation of stellar velocity ellipsoids within the discs of spiral galaxies. The BEADS-2D code directly evolves observable quantities, and it can follow their evolution over many e-folding times starting from arbitrary small perturbations. These properties make the stellar hydrodynamics approach advantageous over the more common N-body codes which, however, suffer strongly from particle noise. Presently, the BEADS-2D code is limited by its two-dimensionality (the thin-disc approximation) but a fully three-dimensional implementation is currently under development. The interested reader is referred to VT06 for details. In the current simulations, the numerical resolution has grid zones that are equally spaced in the -direction and logarithmically spaced in the -direction. The inner and outer reflecting boundaries are at 0.2 kpc and 35 kpc, respectively. We extend the outer computational boundary far enough to ensure a good radial resolution in the region of interest within the inner 20 kpc. For instance, the radial resolution at 1 kpc and 20 kpc is approximately 10 pc and 150 pc, respectively.

## 3 Model galaxy and initial conditions

In this paper we study the influence of a spiral gravitational field on the properties of stellar velocity ellipsoids within the disc of a spiral galaxy. Our model galaxy consists of a thin, self-gravitating stellar disc embedded in a static dark matter halo. The initial surface density of stars is axisymmetric and is distributed exponentially according to

(2) |

with a radial scale length of 4 kpc. The central surface density is set to pc. We note that the densities of stellar discs are indeed found to decay exponentially with distance, with a characteristic scale length increasing from kpc for the early type galaxies to kpc in the late type galaxies (Freeman freeman ).

The gravitational potential of the stellar disc is calculated by

(3) | |||||

This sum is calculated using a FFT technique that applies the 2D Fourier convolution theorem for polar coordinates (Binney & Tremaine BT , Sect. 2.8).

The initial mean rotational (tangential) velocity of stars in the disc is chosen according to

(4) |

The transition radius between the inner region of rigid rotation and a flat rotation in the outer part is given by , which we set to 3 kpc. The smoothness of the transition is controlled by the parameter , set to 3. The velocity at infinity, , is set to 208 km s. The resulting rotation curve is shown in Fig. 1 by the solid line. The initial mean radial velocity of the stars is set to zero.

The radial component of the velocity dispersion is obtained from the relation for a given value of the Toomre parameter . Here, is the epicycle frequency. We assume that throughout most of the disc is constant and equal 1.3 but is steeply increasing with radius at kpc. The initial radial distribution of is shown in Fig. 1 by the dashed line. The tangential component of the velocity dispersion is determined adopting the epicycle approximation, in which the following relation between and holds (Binney & Tremaine BT , p. 125):

(5) |

where is the circular speed, i.e. the speed of a hypothetical star in a circular orbit determined exclusively by the gravitational potential. The value of is yet undefined (since the gravitational potential of the dark matter halo is yet undefined) and we substitute the circular speed with the initial rotation velocity of stars . We emphasize that equation (5) is valid only for stars on nearly circular orbits. We discuss the practical application of the epicycle approximation to spiral galaxies in section 5.1. The mixed velocity dispersion is initially set to zero.

Once the rotation curve and the radial profiles of the stellar surface density and velocity dispersions are fixed, the dark matter halo potential can be derived from the following steady-state momentum equation (see eq. (4) in VT06)

(6) |

## 4 Development of a spiral structure

When starting our simulations we add a small random perturbation to the initially axisymmetric surface density distribution of stars. The relative amplitude of the perturbation in each computational zone is less than or equal to . Since the model stellar disc is characterized by , it becomes vigorously unstable to the growth of non-axisymmetric gravitational instabilities. We find that the positive density perturbation relative to the initial density distribution in the disc can be a good tracer for the spiral structure. Hence, we define the relative stellar density perturbation as

(7) |

where is the initial axisymmetric stellar density distribution in the disc and is the current stellar density distribution. The positive values of at Gyr are plotted in Fig. 2 by the red contour lines. The maximum and minimum contour levels correspond to the relative perturbations of and , respectively, and the step is . The position of corotation at approximately 8.5 kpc is shown by the dashed circle. The stellar disc has clearly developed a two-armed spiral structure and a central bar. A movie showing the time development of the spiral structure can be found at http://www.univie.ac.at/theis/swing_VT07.avi. Snapshots from different stages of the evolution have also been presented by VT06 for a similar model stellar disc.

In our numerical simulations, the bar forms within the 4:1 ultraharmonic resonance (UHR), which is located at approximately 4.5 kpc. At the end of the bar, a slightly weaker spiral structure emerges almost perpendicularly to the bar. Such transitions of bars to spirals are not uncommon (e.g. NGC 1672 or NGC 3513). Radially, our spiral arms reach up to the OLR ( 15 kpc). Azimuthally, they cover an angle of about where they start to fade out. The strongest gradients in the stellar density are found near the outer (convex) edges of the spiral arms. The pattern speeds of the bar and the spiral structure are identical. Based on orbital analysis, Kaufmann & Contopoulos (kaufmann96 ) identified the UHR as a transition region between a bar and spiral features for systems characterized by a single pattern speed for both the bar and the spiral. Recently, Patsis (patsis06 ) has shown that such a coupling is related to chaotic orbits that resemble for a long time the 4:1-resonance orbital behaviour. We note that beyound 10 kpc our spiral arms are very tightly wound due to the decreasing pitch angle. This is in contrast to logarithmic spirals observed in most spiral galaxies. This inconsistency may be caused by the presence of the OLR in our model galaxy where spiral arms tend to converge to a resonant ring. We also emphasize that the amplitude of the spiral stellar density perturbation is quickly decreasing as one approaches the OLR, which makes it difficult to observe the stellar spirals near the OLR in real galaxies.

The non-axisymmetric gravitational
field of the spiral arms and the bar perturbs stellar orbits. To better
illustrate this effect, we plot in Fig. 2 the residual velocity field of stars
at Gyr. The residual velocities are obtained by subtracting the
initial mean rotation velocities of stars (shown in Fig. 1) from the
current mean velocities. Noticeable large-scale perturbations to the initially
axisymmetric flow of stars are seen between the spiral arms and near the
bar. The flow pattern near the convex edges of spiral arms (and outside
corotation) is particularly interesting. The non-axisymmetric gravitational
field of the spiral arm creates two counterpointed large-scale streams
of stars, which merge near the convex edge of the spiral into a narrow
retrograde stellar stream^{3}^{3}3We use the notion “stellar stream” to
describe a large-scale non-circular motion of stars in our numerical simulations.
These streams may not be equivalent to the observed moving groups of stars,
which are characterized by velocities and ages that are different from
their environment. In principle, multicomponent stellar-hydrodynamics or N-body simulations are needed
to compare our stellar streams with the observed moving groups..
A similar phenomenon, though less apparent,
is seen near the inner edges of the bar where two merging streams of stars
produce a strongly perturbed velocity pattern. The regions with merging stellar streams
occupy only a small portion of the stellar disc.
However, as we will see in the next section, these regions are characterised
by the most pronounced influence of the non-axisymmetric gravitational field on the shape
and orientation of stellar velocity ellipsoids.

As demonstrated recently in VT06, the most likely physical interpretation for the growth of a spiral structure in our model disc is swing amplification. Amplification occurs when any leading spiral disturbance unwinds into a trailing one due to differential rotation (Goldreich & Lyndell-Bell GLB , Toomre Toomre , Athanassoula athanassoula84 , and others). According to Julian & Toomre (JT ), the gain of the swing amplifier in a stellar disc with and a flat rotation curve is the largest when . The latter quantity is defined as , where is the circumferential wavelength of an -armed spiral disturbance and is the longest unstable wavelength in a cold disc. Generally speaking, swing amplification is most efficient when . At and , the swing amplifier is strongly moderated by local shear and random motion of stars, respectively.

The top panel of Fig. 3 shows the radial profiles of the -parameter for the (solid line), (dashed line), (dash-dotted line), and (dash-dot-dotted line) spiral perturbations. The region of maximum swing amplification is bound by the two horizontal dotted lines. It can be clearly seen that the perturbations are swing amplified in the inner 5 kpc, whereas the perturbations are swing amplified over a much larger area between 1 kpc and 12 kpc. Consequently, the perturbations are expected to engage more stellar mass into the amplification process than the perturbations. Hence, it is natural to conclude that the density perturbations should grow faster and ultimately surpass the perturbations.

Why then does our model stellar disc develop a two-armed global spiral pattern rather than a three-armed or any other multi-armed pattern? Indeed, higher order perturbations are also swing amplified over a considerable area of the disc, as Fig. 3 (top panel) demonstrates. The answer is that swing amplification needs a feedback mechanism that can constantly feed a disc with leading spiral disturbances. The trailing short-wavelength disturbances may propagate through the centre and emerge on the other side as leading ones, thus providing a feedback for the swing amplifier. Figure 4 in VT06 and the mentioned movie nicely illustrate this phenomenon. However, the propagation of stellar disturbances through the disc centre is possible if there is no inner Lindblad resonance (ILR) or if the waves are reflected before reaching the ILR. A -barrier allows for such a reflection of propagating waves (e.g. Athanassoula athanassoula84 ) but it is absent in the inner region of our model disc (see Fig. 1).

To determine the position of Lindblad resonances , we plot in the bottom panel of Fig. 3 the radial profiles of and for the initial axisymmetric stellar disc and two global spiral modes and . The angular velocity of the two-armed global spiral pattern (see next section) is determined from both visual and numerical estimates to be km s kpc. This value is plotted by the dotted line in the bottom panel of Fig. 3. It is clearly seen that the disturbances have no inner Lindblad resonance, at all radii. On the other hand, the disturbances have the inner Lindblad resonance near the disc centre, in the inner 2 kpc. This prevents the (and any higher order) disturbances from propagating through the disc centre, which terminates the feedback mechanism for the swing amplifier for these modes. We find that this mechanism works nicely in stellar discs with other radial configurations of the -parameter and pattern angular velocities and can allow for the dominant growth of , , and other higher order global spiral modes. The results of this study will be presented in a future paper.

## 5 The shape and orientation of stellar velocity ellipsoids

### 5.1 Stellar velocity dispersions and spiral structure

Despite more exotic means of galactic disc heating there are currently two known mechanisms involving spirals as the heating agents. These are transient spiral stellar density waves (e.g. Wielen Wielen , Carlberg & Sellwood CS , Jenkins & Binney JB , De Simone et al. DeSimone ) and multiple spiral stellar density waves (Minchev & Quillen Minchev06 ). The latter mechanism involves two sets of steady state spiral structure moving at different pattern speeds. We note that a single steady state spiral structure as defined in e.g. Lin et al. (Lin ) does not heat stellar discs.

In this section, we focus on the effect that the non-axisymmetric gravitational field of a swing amplified spiral structure may have on the spatial distribution of stellar velocity dispersions within a stellar disc. Figures 4a and 4b show the spatial distribution of the radial and tangential stellar velocity dispersions, respectively, at Gyr. A two-armed spiral structure and a bar are apparent in the spatial distribution of and . To better illustrate this phenomenon, we plot in Figs 4c and 4d the stellar velocity dispersion perturbations relative to their initial values at Gyr. More specifically, Fig. 4c shows the relative perturbations in the radial velocity dispersion calculated as , where is the initial radial velocity dispersion and is the model’s known radial velocity dispersion at Gyr. To avoid confusions, we note that and are measured in km s throughout the paper. Figure 4d shows the relative perturbations in the tangential velocity dispersion calculated in a similar manner. The positive and negative perturbations are plotted with the shadows of red and blue, respectively. The contour lines of the positive stellar density perturbations are shown for convenience. Figure 4 clearly demonstrates that the non-axisymmetric gravitational field has a pronounced global influence on the spatial distribution of stellar velocity dispersions in the disc. The most distinctive feature in Fig. 4 is that the positive (and negative) relative perturbations in density and both stellar velocity dispersions correlate spatially, suggesting a causal link. We emphasize that the largest spatial gradients in both velocity dispersions (and in the stellar surface density) are found near the outer (convex) edges of the spiral arms. As we will see in Section 6.1, this may cause a severe failure of the epicycle approximation in these regions.

### 5.2 The shape of stellar velocity ellipsoids

In contrast to gaseous discs, stellar discs are essentially non-isotropic. This means that the local properties of stars, such as mean velocities, are determined by the velocity dispersion tensor rather than by an isotropic pressure. The tensor is often non-diagonal in the local coordinate system (). The principal axes of a diagonalized velocity dispersion tensor form an imaginary ellipsoidal surface that is called the velocity ellipsoid. The available measurements in the solar vicinity indicate that the ratio of the smallest versus largest principal axes of the stellar velocity ellipsoid in the disc plane does not vary significantly with the colour. According to Dehnen & Binney (Dehnen ), the axis ratio is approximately for stars with . To avoid confusions, we note that throughout the paper and are measured in km s. The colour of an ensemble of stars may be considered as a rough estimate of their mean age. A weak sensitivity of the axis ratio to the age of stellar populations is intriguing since the principal axes of velocity ellipsoids (when considered separately) are indeed age sensitive and show a noticeable (a factor of two) increase along increasing colours (Dehnen & Binney Dehnen ).

If the axis ratio is indeed only weakly sensitive to the age of stellar populations, it would be interesting to determine other phenomena that could possibly affect the values of . In this section, we consider the effect that the non-axisymmetric gravitational field of the spiral arms and the bar may have on the shape of stellar velocity ellipsoids within the disc of our model galaxy. Figure 5 shows the contour plot of the positive stellar density perturbation (red lines) superimposed on the grey-scale map of in the disc at Gyr. It is evident that the central 5 kpc are characterized by . Most of the intermediate and outer disc regions have , which is comparable to the solar neighbourhood value of 0.64 (Dehnen & Binney Dehnen ). However, parts of the stellar disc that are immediately adjacent to the outer (convex) edges of spiral arms are characterized by a noticeably smaller value of . To better quantify this phenomenon, we calculate the axis ratio along a radial cut at the azimuthal angle , which is shown in Fig. 5 by the dashed line. The resulting ratio is plotted in Fig. 6 by the solid line, whereas the dashed line shows the relative stellar density perturbation obtained along the same radial cut. The area below the dashed line is filled with gray to better identify the position of spiral arms. The dotted line represents the solar neighbourhood value . It is evident that the axis ratio has deep minima near the outer (convex) edges of the spiral arms. The smallest value of is attained at kpc.

The cause for such small values of near the outer edges of spiral arms is unclear. One of the principal velocity dispersions or becomes much smaller than the other, if the mixed velocity dispersion becomes comparable to the squared radial () and tangential () velocity dispersions. We will see in the next section that the outer (convex) edges of the spiral arms are indeed characterized by large values of and associated vertex deviations. These large values of the mixed velocity dispersion are most likely caused by the two large-scale (counterpointed) streams of stars that merge near the outer edge of the spiral arms into a retrograde narrow stellar stream (see Fig. 2). It is not clear that such streams are a general property of spiral galaxies and high-resolution observation are necessary to confirm their existence.

The (possible) existence of regions with abnormally low ratios of can be potentially used to determine the position of stellar spiral arms. This task is not trivial since the spiral arms in external galaxies are usually manifested by strong emission from current star formation, which traces the gaseous response to the underlying stellar spiral arms. There could exist a substantial phase shift between the gaseous and stellar spiral arms. Figure 5 indicates that in the inter-arm region , while near the convex edge of the arm can become as low as 0.25. Such a large difference can be detected. Moreover, our numerical simulations show that a larger contrast is expected for gravitationally more unstable discs. For instance, in a stellar disc the ratio near the convex edges of spiral arms can become as small as but it stays near in the inter-arm region. On the other hand, multi-armed spirals usually show a much smaller contrast in the ratio, possibly due to smaller perturbations to the velocity field of stars. High-resolution measurements of the stellar velocity ellipsoids in grand-design two-armed spiral galaxies are needed to confirm the existence of abnormally low values of the axis ratio .

Finally, it should be stressed that our models are genuine 2D-models. From a higher-order moment analysis Cuddeford & Binney (cuddeford94 ) showed that neglecting the -motion can lead to a systematic overestimate of and, thus, an overestimate of the Oort ratio . Though such 3D-effects will definitely have an impact on the exact Oort ratio, the changes are about 20-25% (cf. e.g. fig. 2 in Cuddeford & Binney). Therefore, we do not expect them to alter our results qualitatively, i.e. the large spatial variations of found in our analysis will be kept. In fact, if is smaller than , the ratio should be overestimated in our 2D-models and the mentioned variation should be even larger.

### 5.3 Vertex deviation

It has long been known from observations that the mixed velocity dispersion within the disc plane is non-zero for stars of many stellar types (see e.g. Dehnen & Binney Dehnen ). The BEADS-2D code directly evolves the mixed velocity dispersion and is particularly suitable for studying the evolution of in galactic discs. Figure 7 shows the map of at Gyr. Positive and negative values of are shown with the shadows of red and of blue, respectively. We emphasize the fact that may be both negative and positive. Only the values between and are plotted to emphasize in the outer disc. It is evident that takes large values near the convex edges of spiral arms and in the inner central region. For instance, can become as large as near the convex edges of spiral arms and near the disc centre. The presence of substantially non-zero implies that the axis ratio may become quite different from the ratio of the tangential versus radial velocity dispersions.

The non-zero indicate that the stellar velocity ellipsoids are not aligned along the local coordinate vector (galactic centre-anticentre direction). The degree of this misalignment is often quantified in terms of a vertex deviation. A classic definition for the vertex deviation is (Binney & Merrifield BM , p. 630)

(8) |

According to this definition, the vertex deviation takes the values between and . In practice, however, the degree of misalignment of the major axis of stellar velocity ellipsoids with the centre-anticentre direction may exceed these bounds, especially in the regions where the epicycle approximation breaks down. In an extreme case, one might consider a velocity distribution with a vanishing non-diagonal element of and exceeding . Then the major principal axis of the velocity ellipsoid will be aligned with the azimuthal direction, yielding . Therefore, we have extended the classic definition (8) as in VT06

According to this re-definition, the vertex deviation can take the values between and . In this paper we extend our previous analysis (VT06) and search for a correlation between the magnitude (and sign) of the vertex deviation and the characteristics of the spiral gravitational field.

Global Fourier amplitudes are often used to quantify the growth rate of spiral perturbations in gravitationally unstable stellar and gaseous discs. These amplitudes are defined as

(9) |

where is the spiral mode, is the total disc mass, and and are the disc inner and outer radii, respectively. The global Fourier amplitudes provide a rough estimate on the total relative amplitude of a spiral density perturbation with mode .

The solid line in Fig. 8a shows the temporal evolution of the dominant global Fourier amplitude (in log units) in our model stellar disc. The amplitude indicates a steady growth of a two-armed spiral structure in the disc from 0.5 Gyr till approximately 1.6 Gyr, when it saturates at . The subsequent evolution of the stellar disc is characterized by a near-constant value of for approximately another 300 Myr, during which the spiral structure is being destroyed by the gravitational torques from spiral arms (this phenomenon is discussed in detail in VT06). It is interesting to see the temporal behaviour of the vertex deviation as the spiral structure grows, saturates, and disperses. The dashed line in Fig. 8a shows the temporal evolution of the mass-weighted vertex deviation (in log units) defined as

(10) |

where is the stellar mass in the -th computational zone, is the vertex deviation in the -th zone, and summation is performed over all computational zones. We plot the mass-weighted vertex deviation rather than its arithmetic average to deemphasize high vertex deviations in regions of low stellar density. A correlation between and is clearly evident – both grow and saturate in a similar manner. The saturated value of the mass-weighted vertex deviation is approximately . To emphasize the correlation between and , we plot in Fig. 8b the mass-weighted vertex deviation versus the global Fourier amplitude (filled circles), both in log units. The least-square fit shown by the solid lines yields a near-linear relation, .

While the mass-weighted vertex deviation appears to correlate globally with
the amplitude of spiral density perturbations,
it is not clear that this correlation should hold locally as well.
Indeed, Fig. 9 shows the vertex deviations in the disc at Gyr.
The positive and negative vertex deviations are plotted with the shadows of red and
blue, respectively. White space corresponds to a near-zero vertex deviation.
The black contour lines in Fig. 9 delineate the spiral pattern by showing
the regions with a positive relative stellar density perturbation .
The minimum and maximum contour
levels correspond to and , respectively.
It is clearly seen that the value (and sign) of the vertex deviation
does not correlate with the amplitude of the positive stellar density perturbation.
In fact, regions with maximum amplitude are characterized by
near-zero or small negative vertex deviations ().
At the same time,
the convex edges of spiral arms, where the stellar density perturbations are small
but the gradients of the density perturbations are large, are characterized by
quite large negative vertex deviations reaching
in the outer disc. Considerable positive vertex deviations up to
are seen between the spiral arms, where the stellar density perturbations are negative.
Somewhat surprisingly, the bar has a near-zero vertex deviation but
the central regions on both sides of the bar are characterized by
the largest vertex deviations up to .
These large values of are off the scale in Fig. 9
(and subsequent figures) to emphasize the regions with lower vertex deviations^{4}^{4}4We note that the sign of the vertex deviation changes when the galactic rotation is reversed
from the counterclockwise rotation (as in present numerical simulations) to the
clockwise rotation. This happens because the mean tangential velocities ()
and instantaneous tangential velocities () of stars change the sign but
the mean and instantaneous radial velocities ( and , respectively)
do not. As a consequence, the mixed velocity dispersion
defined by equation (1) and the corresponding vertex deviations
change the sign, too..

## 6 Discussion

### 6.1 Comparison with the epicycle approximation

It is often observationally difficult to measure both radial and tangential velocity dispersions of stars. Hence, it is tempting to use equation (5) to derive either of the two dispersions from the known one and the circular speed of stars. One has to keep in mind that equation (5) is valid only in the epicycle approximation, when the epicycle amplitude (deviation from a purely circular orbit) at a specific radial distance is much less than both and the scale length for changes in both, the stellar density and dispersion (Kuijken & Tremaine KT , Dehnen dehnen99 ). The latter condition may become quite stringent, especially near the convex edges of spiral arms. Therefore, it is important to know if equation (5) holds in spiral galaxies with a well-developed spiral structure.

To test the validity of equation (5), we calculate the ratio from the model’s known velocity dispersions and compare it with the ratio derived using equation (5). From the observational point of view, it is usually difficult to know the exact value of the circular speed that enters equation (5), since it requires independent measurements of the axisymmetric gravitational potential in a galaxy. The observers usually apply the local and azimuthally averaged rotation velocities of stars and , respectively, as a proxy for . The validity of this assumption is checked below.

To quantify the deviations from the epicycle approximation in our model spiral galaxy, we calculate the relative errors

(11) |

using different approximations for the circular speed . In particular, Fig. 10 and Fig. 11 show the relative errors that are obtained assuming the local rotational velocity and azimuthally averaged rotational velocity as a substitute for the circular speed , respectively. The positive and negative errors are plotted with the shadows of blue and red, respectively. The spiral pattern in the disc is outlined by contour lines showing the positive relative stellar density perturbations at Gyr.

A visual inspection of Fig. 10 and Fig. 11 indicates that the deviations from the epicycle approximation are largest near the convex edges of spiral arms. To better illustrate this effect, we take a radial cut at in Fig. 10 and Fig. 11 (shown by the dashed line) and plot the resulting radial distribution of by the solid lines in Fig. 12a and Fig. 12b, respectively. The positive relative stellar density perturbation along the radial cut is shown by the dashed lines. The area below these lines is filled with grey to identify the position of spiral arms. It is evident that if the local rotation velocity of stars is used as a proxy for the circular speed (Figure 12a), then the relative errors can become as large as and near the outer (convex) edges of spiral arms. If the azimuthally averaged rotational velocity of stars is used as a proxy for (Figure 12b), the corresponding relative errors become roughly a factor of three smaller. Yet the relative error can become as large as near the outer edge of the innermost spiral at kpc, which renders the epicycle approximation inapplicable in these regions.

Finally, we check the validity of the epicycle approximation when the circular speed is determined from the axisymmetric gravitational potential using equation (6) with and set to zero (cold stellar disc). This approach is accurate but difficult to realize observationally for the reasons explained above. Nevertheless, it might be instructive to compare the ratios of tangential versus radial velocity dispersions obtained in our numerical model with those obtained from the exact epicycle approximation. The corresponding relative errors are shown in Fig. 13 with the shadows of blue (positive errors) and red (negative errors). The spiral pattern in the disc is outlined by the contour lines showing the positive relative stellar density perturbations at Gyr. The convex edges of spiral arms are again characterized by large relative errors. This effect is illustrated in Fig. 12c, which shows (solid line) and (dashed line) calculated along the radial cut at (Fig. 13, dashed line). The relative errors can become as large as near the convex edge of the innermost spiral at kpc and are below elsewhere.

We conclude that the epicycle approximation is severely in error near the convex edges of spiral arms. On the other hand, the deviation from the epicycle approximation in the inner disc ( kpc) remains reasonably small, within per cent. We argue that one should use the azimuthally averaged mean tangential velocity of stars (if available) rather than the local mean tangential velocity as a proxy for and use the epicycle approximation with an extreme caution near the spiral arms. It is interesting to note that the epicycle approximation holds, to within a 20 per cent error, in and near the bar (except for two small spots on both sides of the bar where the relative error can become as large as ).

### 6.2 Origin of vertex deviation

It is generally thought that two factors contribute to the vertex deviation: (i) moving groups and (ii) a large-scale non-axisymmetric component of the galactic gravitational field (Binney & Merrifield BM ). However, it is not clear which of the two factors is more important. The situation may become even more complicated, since the moving groups may partly be caused by the non-axisymmetry in the galactic gravitational field. In order to get a better understanding of the origin of the vertex deviation in our simulation, we analyse the correlation of the vertex deviation with different dynamical quantities.

#### 6.2.1 Correlation of vertex deviation with the gravitational field

We start by examining the role of the non-axisymmetric gravitational field in generating a vertex deviation. Kuijken & Tremaine (KT ) considered non-axisymmetric gravitational potentials typical for a triaxial (ellipsoidal) dark matter halo. They found that the magnitude of the vertex deviation should be proportional to the product of and the ratio of the non-axisymmetric versus axisymmetric contributions to the total gravitational potential . Here is the position angle with respect to the minor axis of the equipotential surfaces of an ellipsoidal halo (see equations [14] and [15] in Kuijken & Tremaine KT ). The vertex deviation is expected to be zero near the minor axis of an ellipsoidal halo . Of course, a direct comparison of our numerical simulations with the analytical predictions of Kuijken & Tremaine is problematic, since the dark matter halo is axisymmetric in our numerical simulations and the non-axisymmetry is introduced by the non-stationary stellar spiral arms and the bar. Moreover, the simulations allow for a completely non-linear evolution, whereas Kuijken & Tremaine’s analysis is based on linear perturbation theory. Nevertheless, it may be helpful to compare two different non-axisymmetric structures and their effect on the spatial distribution of vertex deviations in a stellar disc.

In our model, we calculate the non-axisymmetric gravitational potential by subtracting the gravitational potential of the initial axisymmetric stellar disc from the current model’s known gravitational potential (which is the sum of contributions from the axisymmetric and non-axisymmetric parts of the stellar disc). The ratio of the non-axisymmetric versus axisymmetric contributions to the total gravitational potential in our model is then determined as

(12) |

This quantity calculated at Gyr is shown in Fig. 14 by the contour lines. The minimum and maximum contour levels are and , respectively, and the step is 0.02. We use dashed and solid lines to differentiate between the negative and positive values of , respectively. The positive and negative vertex deviations are shown with the shadows of red and blue, respectively. It is obvious that the magnitude of the vertex deviations in our model stellar disc does not correlate spatially with the ratio of the non-axisymmetric versus axisymmetric stellar gravitational potentials. The maximal vertex deviations (by absolute value) in our model are found along the minor axis of the central bar, whereas in the case of an elliptical halo Kuijken & Tremaine (KT ) predicted near-zero vertex deviations there. The reason for this obvious disagreement is not clear and numerical simulations of vertex deviations in models with static elliptical dark matter halos will be presented in a future paper.

However, we should keep in mind that the dynamics of a star is controlled by the gradient of the gravitational potential rather than by the potential itself. We calculate the absolute value of the stellar gravitational potential gradient in the disc as

(13) |

This quantity can be regarded as the “strength” of the stellar gravitational field in the disc. We further denote the absolute value of the initial stellar gravitational potential gradient and calculate the relative perturbation in the stellar gravitational potential gradient as

(14) |

The spatial distribution of at Gyr is shown in Figure 1 with the contour lines. In particular, the solid lines delineate positive relative perturbations in the stellar gravitational potential gradient and dashed lines delineate negative relative perturbations. The minimum and maximum contour levels are and , respectively. The positive and negative vertex deviations are shown with the shadows of red and blue, respectively. A visual inspection of Fig. 15 reveals a strong correlation between and positive in the inner disc and a mild correlation between these quantities in the outer disc.

To better illustrate this correlation, we find along a radial cut at (dashed line in Fig. 15) and plot the resulting values in Fig. 16 with the dashed line. The corresponding vertex deviation is plotted with the solid line. We should keep in mind that changes its sign when the galactic rotation is reversed but does not. Hence, we show the absolute (rather than actual) value of the vertex deviation. The radial position of the spiral arms and the bar at is outlined by the dotted line which represents the positive relative perturbation in the stellar surface density (the area below this line is filled with grey). The strongest correlation between and is seen in the central region where both quantities attain their maximum values. A noticeable correlation between and is also seen near the outer (convex) edges of spiral arms at kpc and kpc. On the other hand, there is little or no correlation in the inter-arm region at kpc, where is negligible but can become as large as . A similar lack of correlation is seen at kpc. We conclude that the relative perturbation in the gravitational potential gradient can only partly account for the magnitude of the vertex deviation in our model disc.

#### 6.2.2 Correlation of vertex deviation with the streaming motions

In this section, we search for the correlation between the vertex deviation and local properties of the stellar velocity field, in particular, the non-circular streaming motions of stars clearly seen in Fig. 2. The streaming motion of stars can be quantified by calculating the spatial derivatives of the radial and tangential mean stellar velocities, namely, , , and . We have examined all three characteristics of the streaming motion and found that is best suited for our purposes. Figure 17a shows the spatial distribution of (in units of km s kpc) in our model stellar disc at Gyr. The positive and negative values of this quantity are plotted with the shadows of red and blue, respectively. The spatial distribution of the vertex deviation at Gyr, initially shown in Fig. 9, is reproduced in Fig. 17b to facilitate the comparison between and . A striking similarity between the spatial distributions of and is clearly seen in the figure, especially in the outer parts of the stellar disc.

To better illustrate the spatial correlation between and , we take two radial cuts at and (shown with dashed lines in Fig. 17) and plot the resulting radial profiles of (dashed lines) and (solid lines) in Fig. 18. Now both quantities happen to have the same sign and the actual values of the vertex deviation are shown. The comparison of Fig. 16 and Fig. 18 demonstrates that the correlation between and is in general better than the correlation between and . Yet there are regions where and show no correlation. For instance, the derivative is maximal (by negative value) at approximately 14.5 kpc (Fig. 18a) but the corresponding vertex deviation changes its sign and is characterized by negligibly small values there. Other characteristics of the streaming motion ( and ) show similar correlation with the vertex deviation, though less pronounced than in the above case with .

To summarize, the vertex deviation is better correlated with the characteristics
of stellar streams such as than with
the relative perturbation in the gravitational potential gradient in the stellar
disc^{5}^{5}5We stress that the vertex
deviation changes its sign when the rotation is reversed but
the spatial derivatives of the mean stellar velocities do not. Hence, for a clockwise rotation,
we would have obtained an anti-correlation between
and ..
The relative perturbation in the gravitational potential shows essentially no correlation with the vertex deviation, which could be due to the fact
that the dynamics of stars is controlled by the gradient of the gravitational
potential rather than the gravitational potential itself.
It appears that the non-axisymmetric gravitational field of the spiral arms and the bar perturbs stellar
orbits by creating large-scale stellar streams. These streams introduce (under specific
yet poorly understood conditions) the mixed velocity dispersion into
the stellar velocity dispersion tensor and consequently generate the vertex deviation.

This result is in nice qualitative agreement with the linear perturbation analysis in KT91. For instance, their equation (19) shows also a strong correlation between the vertex deviation and the mean velocity field (and its gradients). Moreover, the inner bar region is characterized by rigid rotation, which means KT91’s becomes unity and the denominator formally vanishes yielding a divergent vertex deviation. Such large vertex deviations are in good agreement with the highest values in our simulations found near the central bar.

#### 6.2.3 Higher-order moments

When applying the Boltzmann moment equations each equation of -th order includes a term of -order, by this formally creating an infinite set of equations. In order to end up with a finite set of equations, one has to introduce a closure condition. In hydrodynamics this can be a proper equation of state (e.g. a polytropic equation of state). In stellar hydrodynamics the situation is more complex, because an obvious equation of state does not exist. For instance, an isotropic pressure (or velocity distribution) is guaranteed in a ”normal” gas by collisional relaxation which is usually a rather fast process compared to other timescales like the dynamical timescale. Contrary to a gaseous system, two-body relaxation in a three-dimensional stellar system is a rather slow process. Adopting Spitzer & Hart’s (spitzer71 ) formula for the relaxation time ,

(15) |

we get with a velocity dispersion of 10 km s, a typical
stellar mass of 1 M, a local mass density
and a number of stars^{6}^{6}6Equation (15)
is derived for a spherical system. Thus we
consider a sphere with a radius equal to a typical disc scale height of 100 pc
containing about stars. Though this is only a very rough estimate, it
does not alter our argument, because only shows up in the Coulomb logarithm
term. of , a relaxation time of about
.
This timescale exceeds the dynamical timescale (and the Hubble time) by far.
Even if we adopt the smallest values of about 7 km s for the
velocity dispersions of young stars, the relaxation time is still very long.
Thus, in a purely stellar galactic disc two-body relaxation is negligible,
and heat conduction terms which are related to the third
order moments can be neglected.

In order to check the third order moments quantitatively we started to perform numerical test particle simulations in a given time-dependent disc potential including a spiral. The numerical treatment is similar to the one adopted in Blitz & Spergel (blitz91 ) and more recently in Minchev & Quillen (minchev07a ). First results of these simulations show that the 3rd order terms are small compared to the other moments, by this confirming our numerical treatment. A more detailed analysis of these simulations will be given elsewhere.

## 7 Summary and Conclusions

In this paper, we present two-dimensional stellar hydrodynamics simulations of a flat galactic disc embedded in a standard dark matter halo, yielding a rigid rotation curve in the central part and a flat rotation curve outside 3 kpc. We solve the Boltzmann moment equations up to second order in the thin-disc approximation and follow the growth of a two-armed spiral structure deeply into the non-linear regime. We analyse in detail the shape and orientation (vertex deviation) of the stellar velocity ellipsoids when the spiral structure has reached a saturation phase. The simulations yield the following main results.

(i) The shape of the stellar velocity ellipsoids (as defined by the ratio of the smallest versus largest axes of the stellar velocity ellipsoid) and vertex deviations show large variations over the entire stellar disc. In particular, the radial distribution of the axis ratio shows deep minima near the outer (convex) edges of stellar spiral arms.

(ii) The spatial distributions of the radial and tangential stellar velocity dispersions show a similar two-armed symmetry as the underlying distribution of the stellar surface density.

(iii) The epicycle approximation fails in the regions affected by stellar spiral arms. In particular, it cannot reproduce the tangential-to-radial stellar velocity dispersion ratio obtained from the numerical simulations. The deviations from the epicycle approximation are especially large near the outer (convex) edges of stellar spiral arms, where the epicycle approximation may be in error by a factor of two or more.

(iv) Globally, there is an excellent correlation between the mass-weighted vertex deviation (integrated over the whole stellar disc) and the global Fourier amplitude of the stellar density perturbations. However, locally, there is no simple correlation between the vertex deviation and the mass distribution.

(v) The magnitude of the vertex deviation is not correlated with the gravitational potential (mass distribution) and is only weakly correlated with the gradient of the gravitational potential. However, the vertex deviation correlates strongly with the spatial gradients of mean stellar velocities, in particular, with the radial gradient of the mean radial velocity,

(vi) Large-scale, non-circular stellar streams induced by the non-axisymmetric gravitational field of the stellar spiral arms and the bar may be responsible for the failure of the epicycle approximation, large vertex deviations, and peculiar shapes of the stellar velocity ellipsoids.

Our results are in good qualitative agreement with the results of the linear perturbation analysis given in KT91. Moreover, our simulations include also the non-linear stages of the disc evolution when dynamical feedback limits the growth of a spiral structure. The large spatial variations of the properties of the stellar velocity distribution found in our numerical simulations and the mismatch to the epicycle approximation are a caveat for any analysis of observational data which cover larger spatial areas or which include stellar samples coming from a larger spatial region. It will be interesting to check whether features like the ”u-anomaly” could also be created by large-scale streams induced by a transient or swing-amplified spiral.

## Acknowledgements

We thank the referee, Ivan Minchev, for useful suggestions that helped improve the manuscript. This project was partly supported by grants RFBR 06-02-16819-a, South Federal University grant 05/6-30, and Federal Agency of Education (project code RNP 2.1.1.3483). E.I.V. acknowledges support from a CITA National Fellowship. Numerical simulations were done on the Shared Hierarchical Academic Research Computing Network (SHARCNET). The authors are grateful to Panos Patsis for stimulating discussions.

## References

- 2005 Alcobé S., Cubarsi R., 2005, A&A, 442, 929
- 1984 Athanassoula E., 1984, Phys. Rep., 114, 319
- 2000 Bertin G., 2000, Dynamics of Galaxies, Cambridge Univ. Press
- 1999 Bienaymé O., 1999, A&A, 341, 86
- 1998 Binney J., Merrifield M. S., 1998, Galactic Astronomy, Princeton Univ. Press
- 1987 Binney J., Tremaine S., 1987, Galactic Dynamics, Princeton Univ. Press
- 1991 Blitz L., Spergel D.N., 1991, ApJ, 370, 205
- 1985 Carlberg R.G., Sellwood J.A., 1985, ApJ, 292, 79
- 1994 Cuddeford P., Binney J., 1994, MNRAS, 266, 273
- 1998 Dehnen W., 1998, AJ, 115, 2384
- 1999 Dehnen W., 1999, AJ, 118, 1190
- 1999 Dehnen W., 2000, ApJ, 119, 800
- 1998 Dehnen W., Binney J.J., 1998, MNRAS, 298, 387
- 2003 Mühlbauer G., Dehnen W., 2003, A&A, 401, 975
- 2004 De Simone, R., Wu, X., Tremaine, S., 2004 MNRAS, 350, 627
- 1906 Eddington A.S., 1906, MNRAS, 67, 34
- 1997 ESA, 1997, The Hipparcos and Tycho Catalogues (ESA Sp-1200) (Nordwijk ESA)
- 1970 Freeman K.C., 1970, ApJ, 160, 811
- 1997 Gerssen J., Kuijken K., Merrifield M.P., 1997, MNRAS, 288, 618
- 1965 Goldreich, P, Lynden-Bell, D, 1965, MNRAS, 130, 125
- 2005 Hogg D.W., Blanton M.R., Roweis S.T., Johnston K.V., ApJ, 629, 268
- 1990 Jenkins, A., Binney, J., 1990, MNRAS, 245, 305
- 1966 Julian W.H., Toomre A., 1966, ApJ, 146, 810
- 1905 Kapteyn J.C., 1905, Report of the Brit. Ass. for the Advancement of Science, p. 257
- 1996 Kaufmann D.E., Contopoulos G., 1996, A&A, 309, 381
- 1986 Kerr F.J., Lynden-Bell D., 1985, MNRAS, 221, 1023
- 1890 Kobold H., 1890, Astron. Nachr. 125, 65
- 1991 Kuijken K., Tremaine S., 1991, in Proc. of Dynamics of Disc Galaxies, Göteborg, B. Sundelius (ed.), p. 71 (KT91)
- 1994 Kuijken K., Tremaine S., 1994, ApJ, 421, 178
- 1969 Lin, C. C., Yuan, C., Shu, F. H., 1969, ApJ, 155, 721
- 1970 Mayor M., 1970, A&A, 6, 60
- 1972 Mayor M., 1972, A&A, 18, 97
- 2006 Minchev, I., Quillen, A. C., 2006, MNRAS, 368, 623
- 2007 Minchev I., Quillen A.C., 2007, MNRAS, 377, 1163
- 2007 Minchev I., Nordhaus J., Quillen A.C., 2007, ApJ, 664, L31
- 2006 Patsis P., 2006, MNRAS, 369, L56
- 1997 Polyachenko V. L., Polyachenko E. V., Strel’nikov A. V., 1997, Astron. Zhurnal, 23, 598 (translated Astron. Lett. 23, 525)
- 1997 Ratnatunga K.A., Upgren A.R., 1997, ApJ, 476, 811
- 1907 Schwarzschild K., 1907, Nachr. Kgl. Ges. d. Wiss. zu Göttingen, Math. Phys. Klasse, 5, 614
- 1971 Spitzer L., Hart M.H., 1971, ApJ, 164, 399
- 1981 Toomre, A., 1981, in Fall S. M., Lynden-Bell D., eds, The Structure and Evolution of Normal Galaxies. Cambridge Univ. Press, Cambridge
- 2006 Vorobyov E.I., Theis Ch., 2006, MNRAS, 373, 197 (VT06)
- 1974 Wielen R., 1974, Highlights of Astronomy, 3, 395
- 1977 Wielen R., 1977, A&A, 60, 263
- 1971 Yuan C., 1971, AJ, 76, 664