© ESO 2021

1. Introduction

The ubiquity of grand-design bisymmetric spiral arms that are observed in galactic disks in the K band indicates that the paradigm of the “density wave” remains a valid dynamical model with which the spiral structure of many disk galaxies can be described. While the scenario of quasi-stationary spiral density waves has been strongly criticized (see Toomre 1977; Athanassoula 1984; Binney & Tremaine 2008; Sellwood & Carlberg 1984; Sellwood 2010; Dobbs & Baba 2014 for reviews), the existence of mechanisms regulating the growth rate of spiral density waves, and hence ensuring their longevity for several pattern periods, remains largely an open question (Bertin 1980; Donner & Thomasson 1994, see Bertin & Lin 1996; Contopoulos 2002 for reviews). The density wave theory fits the description of spiral arms betten when the spiral amplitude does not exceed a value of 10%−20% over a timescale of a few (≈5) pattern rotations. In this regime, nonlinear corrections (Contopoulos 1970; Vandervoort 1971; Norman 1978) to the linear Lindblad-Lin-Shu theory are required (Lindblad 1940, 1961; Lin & Shu 1964, 1966).

The basic formulation of density wave theory (see Binney & Tremaine 2008) deals with spiral wave perturbations in an axisymmetric disk that is regarded in the framework of either the Boltzmann equation for collisionless matter or the hydrodynamical equations for gas. On the other hand, Lindblad (1955) pioneered the orbital description of spiral density waves. In this description, the central objective is to identify the main families of stellar orbits that support the spirals, and also to show how these orbits can collaborate to match the imposed with the response model of spiral arms (see Contopoulos 1971; Berman & Mark 1977; Monet & Vandervoort 1978; Contopoulos & Grosbøl 1986; Patsis et al. 1991). In the orbital description of spiral density waves, a key element are approximately elliptical stable periodic orbits whose orientation changes with the distance from the center (parameterized by the Jacobi energy) so as to produce a response density that closely follows the minima of the imposed spiral potential. The importance of these periodic orbits as the building blocks of spiral arms was first stressed by Lindblad (1956, 1957, 1958, 1960, 1961), who called them “dispersion orbits”. Contopoulos (1970, 1975) developed the theory of resonant periodic orbits near the inner Linblad resonance, providing the main formulas that allow predicting the number and stability of the periodic orbits as a function of the Jacobi energy. The resonant theory of Contopoulos, based on epicyclic action-angle variables, also allows predicting the structure of the phase space around the periodic orbits, yielding the corresponding invariant curves as the level curves of a resonant “third integral” (see Contopoulos 2002 for a review). The predictions of this theory were confirmed by Vandervoort (1973, 1975; see also Vandervoort 1975 and Monet & Vandervoort 1978), Mertzanides (1976) and Berry & Smet (1979). These studies provided figures of the corresponding phase portraits around the elliptical closed orbits computed analytically (by the third integral) or numerically in simple models of galactic potentials with a spiral perturbation. On the other hand, Kalnajs (1973) popularized the orbital version of the density wave theory by providing a schematic figure of the emergence of spiral arms due to the change in orientation of the axes of the elliptical orbits (“precessing ellipses”), a figure that was later reproduced in many reviews and books on the subject. Kalnajs (1973) also computed an approximative formula for the response spiral potential due to the precessing ellipses, which holds beyond the Wentzel-Kramers-Brillouin (WKB) limit of tightly wound spirals. The overall limit of the applicability of the precessing ellipses was considered by Contopoulos (1985). The change in orientation of the elliptical orbits with the energy is explained by resonant perturbation theory (Contopoulos 1975, 2002; Monet & Vandervoort 1978, see a tutorial presentation in Efthymiopoulos 2010). Numerical examples of spirals supported by precessing ellipses have been computed in models of spiral galaxies (Contopoulos & Grosbøl 1986; Patsis et al. 1991; Patsis & Grosbøl 1996; Pichardo et al. 2003; Efthymiopoulos 2010; Tsigaridi & Patsis 2013; Pérez-Villegas et al. 2015; Chaves-Velasquez et al. 2019).

An important outcome of linear density wave theory, which is valid for small amplitudes of the spiral perturbation, is the dispersion relation (see Binney & Tremaine 2008, p. 493), which relates the frequency (pattern speed Ω_sp of the spiral arms) with the radial wavenumber k. For tightly wound m = 2 logarithmic spirals, the radial wavenumber k and the pitch angle α are related by α = cot⁻¹(k r/2). In order to understand the relevance of the dispersion relation to our results below, we point out that in a given (and fixed) axisymmetric model, the dispersion relation depends on the quantities κ(r) (epicyclic frequency), Ω(r) (angular velocity of the circular orbit of radius r), and Σ₀(r) (axisymmetric disk, as well as on the velocity dispersion σ(r) at radius r. All these quantities are completely specified by the model that is used to represent the axisymmetric component of the disk, that is, by the potential V_ax(r) and by the form Q(r) of the profile of Toomre’s Q-parameter as a function of the distance from the center (see Eq. (15) below). Thus, after an axisymmetric model is fixed, the dispersion relation should be regarded as a relation allowing us to specify the form of the function k(r; Ω_sp), where k = dϕ_m(r)/dr and ϕ_m is the phase (position of maxima) of the m-fold spiral mode. Hence, solving the differential equation (for m = 2 spiral arms) dϕ₂(r)/dr = k(r; Ω_sp), we obtain the form of the phase ϕ₂(r; Ω_sp) corresponding to the maxima of the spiral arms. This implies that within a fixed axisymmetric background (potential + velocity dispersion), the linear dispersion relation establishes a correlation between the form of the spiral arms (given by the function ϕ₂) and the value of the spiral pattern speed Ω_sp. It is important to emphasize that in the linear regime, this relation is independent of the amplitude of the spiral perturbation, which is simplified altogether from the equations of the linear theory.

However, in nonlinear wave theory, the amplitude of the spiral perturbation A(r) also has to be taken into account (see Contopoulos 1975, 1980; Norman 1978). This leads to a form of the dispersion relation that correlates all three basic parameters of the spirals, namely, the amplitude A(r), the pitch angle a, and the pattern speed Ω_sp.

Regarded from the orbital point of view, two important nonlinear effects are introduced when the amplitude of the spirals becomes large:

• (i)

  the ‘precessing ellipses’ (closed orbits of the x₁ family) become largely distorted with respect to the elliptical shape, and

• (ii)

  chaos is introduced in the system, leading to the lack of a sufficient number of regular orbits to support the spiral arms.

While these two phenomena have been stressed in previous studies (Contopoulos 1985; Contopoulos & Grosbøl 1986; Patsis et al. 1991, see also Quillen & Minchev 2005), a systematic study yielding the permissible space in all three parameters (A, a, and Ω_sp) within which the precessing ellipses can support spiral arms is still lacking, in the literature. This constitutes the primary motivation for our present study. In particular, we first work with a simple model in which the axisymmetric components (disk and halo) have parameters with values that are pertinent to a Milky Way model. A bulge is added to properly secure good kinematic properties of the orbits at all distances across the disk. We then superpose a m = 2 logarithmic spiral term to the axisymmetric potential with a varying amplitude, pitch angle, and pattern speed, and study the corresponding phase-space structure in the region from the inner Lindblad resonance (ILR) to corotation. Our study focuses on the structure of the characteristic curve and the stability of the x₁ family of periodic orbits, along with the extent of the domain of regular orbits around the x₁ family. We also discuss the role of other families generated near the ILR (in particular, the x₂ family), as well as the range in the parameter values in which the change in pericenter over galactocentric distance r is such as to support a spiral response similar to the one imposed by the model.

Second, stemming from the basic steps of resonant perturbation theory at the ILR, a byproduct of our analysis is to show the convenience of a particular set of variables in illustrating all phenomena relevant to the phase-space structure in the region from the IRL and up to the 4:1 resonance. Thus, we depict all numerical phase portraits (Poincaré surfaces of section) using the epicyclic set of canonical variables (ξ, p_ξ), where ξ = r − r_c, $p_{\xi }=\dot{\xi }$, and r_c is the radius of the circular orbit of a fixed Jacobi energy h in the flow under the axisymmetric potential alone. Furthermore, we depict the characteristic curves of closed orbits (x₁, x₂, etc.) showing the quantity

$$\begin{array}{c}\hfill S(r_c)={({\xi }_0^2+\frac{p_{\xi ,0}^2}{{\kappa }^2(r_c)})}^{1/2}\end{array}$$(1)

against r_c, where (ξ₀, p_ξ, 0) is the fixed point at which a closed orbit intersects a surface of section defined by a fixed azimuth. This choice of variables (instead of the customary choice, i.e., ξ₀ as a function of h) is motivated by the fact that according to resonant perturbation theory, the relevant quantity characterizing the precessing ellipses is the amplitude of the epicyclic oscillation for a particular closed orbit, quantified by S(r_c), as a function of the average distance of the orbit from the center of the disk. We find that spiral arms are in general supported by orbits for which S(r_c) is a decaying function of r_c (see Sect. 4).

Third, we identify the appearance of chaos as the main dynamical phenomenon defining the boundary of the domain in parameter space (A, a, and Ω_sp), outside which spiral arms supported by precessing ellipses cannot exist. A good criterion for determining this boundary is obtained by computing the threshold in (A, a, and Ω_sp) beyond which the x₁ family itself becomes unstable at distances already close to the ILR. This threshold practically marks the nearly complete disappearance of ordered orbits in phase space. Detailed numerical evidence of these phenomena is given at the end of Sect. 4, leading to conclusions that approximately agree with the predictions of density wave theory, and probably also with observations (see the discussion at the end of Sect. 4). In particular, this orbital analysis supports the prediction that larger spiral amplitudes A as well as slower rotation (smaller Ω_sp) are consistent with more open spirals (higher values of the pitch angle).

The paper is structured as follows. Section 2 presents our model. Section 3 discusses the construction of spiral density waves using closed orbits. In Sect. 4, an explanation is given of how periodic orbits are determined and how the phase space is constructed. Moreover, the main results of the paper are presented in this section, where we study the role of the amplitude of the spiral perturbation, the pattern speed of the spirals, and the pitch angle in creating realistic spiral density waves. Finally, we summarize our conclusions in Sect. 5.

2. Model

We considered a model of a spiral galaxy that contains a combination of an axisymmetric and a spiral potential,

$$\begin{array}{c}\hfill V=V_{\mathrm{ax}}+V_{\mathrm{sp}}.\end{array}$$(2)

The axisymmetric potential V_ax is composed of a disk, a bulge, and a halo,

$$\begin{array}{c}\hfill V_{\mathrm{ax}}=V_{\mathrm{d}}+V_{\mathrm{b}}+V_{\mathrm{h}}.\end{array}$$(3)

For the disk potential V_d, we used a Miyamoto-Nagai model (Miyamoto & Nagai 1975) given by the relation

$$\begin{array}{c}\hfill V_{\mathrm{d}}=\frac{-G{M}_{\mathrm{d}}}{\sqrt{r^2+{(a_{\mathrm{d}}+\sqrt{z^2+b_{\mathrm{d}}^2})}^2}},\end{array}$$(4)

where M_d = 8.56 × 10¹⁰ M_⊙ is the total mass of the disk, a_d = 5.3 kpc and b_d = 0.25 kpc. In order to have a 2D disk model, we took z = 0 and $r=\sqrt{x^2+{\mathit{y}}^2}$. For the bulge, we used a Plummer potential V_b given by the relation

$$\begin{array}{c}\hfill V_{\mathrm{b}}=\frac{-G{M}_{\mathrm{b}}}{\sqrt{r^2+b^2}},\end{array}$$(5)

where M_b = 5 × 10¹⁰ M_⊙ is the total mass of the bulge, $r=\sqrt{x^2+{\mathit{y}}^2}$ and b = 1.5 kpc.

The halo potential was a γ-model (Dehnen 1993) with parameters as in Pettitt et al. (2014),

$$\begin{array}{c}\hfill V_{\mathrm{h}}=\frac{-G{M}_{\mathrm{h}(\mathrm{r})}}{r}-\frac{-G{M}_{\mathrm{h},0}}{\gamma r_{\mathrm{h}}}{[-\frac{\gamma }{1+{(r/r_{\mathrm{h}})}^{\gamma }}+ln{(1+\frac{r}{r_{\mathrm{h}}})}^{\gamma }]}_r^{r_{\mathrm{h},\max }},\end{array}$$(6)

where r_h, max = 100 kpc, γ = 1.02, and M_h, 0 = 10.7 × 10¹⁰ M_⊙, and M_h(r) was given by the function:

$$\begin{array}{c}\hfill M_{\mathrm{h}(\mathrm{r})}=\frac{M_{\mathrm{h},0}{(r/r_{\mathrm{h}})}^{\gamma +1}}{1+{(r/r_{\mathrm{h}})}^{\gamma }}.\end{array}$$(7)

The spiral potential is given by the value V_sp for z = 0 of the 3D logarithmic spiral model V_sp(r, ϕ, z) introduced by Cox & Gómez (2002) (see formula (19) in Efthymiopoulos et al. 2020). We have on the disk plane

$$\begin{array}{c}\hfill V_{\mathrm{sp}}=4\pi G{h}_{\mathrm{z}}{\rho }_{0}G(r)exp(-\left(\frac{r-r_0}{R_{\mathrm{s}}}\right))\frac{C}{\mathit{KB}}cos\left[2(\phi -\frac{ln(r/r_0)}{tan(\alpha )})\right],\end{array}$$(8)

where

$$\begin{array}{c}\hfill K=\frac{2}{r|sin(\alpha )|},B=\frac{1+K{h}_z+0.3{(K{h}_z)}^2}{1+0.3K{h}_z},\end{array}$$(9)

and C = 8/(3π), h_z = 0.18 kpc, r₀ = 8 kpc, R_s = 7 kpc, and a = −13° was the pitch angle of the spiral arms. The function G(r) plays the role of a smooth envelope determining the radius beyond which the spiral arms are important. We adopted the form G(r) = b − carctan((R_s0 − r)), with R_s0 = 6 kpc, b = 0.474, and c = 0.335. The spiral density was ρ₀ = 5 × 10⁷, 15 × 10⁷, or 30 × 10⁷ M_⊙ kpc⁻³ in the three different models under study, respectively. These three values of the density were chosen so as to yield spiral F-strength values consistent with those reported in the literature for a weak intermediate and strong spiral, respectively (see, e.g., Block et al. 2004).

The F-strength (Buta et al. 2009) can be defined as either the ratio of the maximum tangential force of the spiral perturbation over the radial force of the axisymmetric background,

$$\begin{array}{c}\hfill F_{\tan }(r)=\frac{{\langle F_{\mathrm{sp}}^{\tan }(r)\rangle }_{{}_{\max }}}{F_{\mathrm{r}}(r)}=\frac{{\langle \frac{1}{r}\frac{\partial V_{\mathrm{sp}}}{\partial \phi }\rangle }_{\max }}{\frac{\partial V_{\mathrm{ax}}}{\partial r}},\end{array}$$(10)

or the ratio of the maximum total force of the spiral perturbation over the radial force of the axisymmetric background, given by the relation

$$\begin{array}{c}\hfill F_{\mathrm{all}}(r)=\frac{\langle F_{\mathrm{sp}}(r)\rangle }{F_r(r)}=\frac{{\langle \sqrt{{\left(\frac{1}{r}\frac{\partial V_{\mathrm{sp}}}{\partial \phi }\right)}^2+{\left(\frac{\partial V_{\mathrm{sp}}}{\partial r}\right)}^2}\rangle }_{\max }}{\frac{\partial V_{\mathrm{ax}}}{\partial r}}.\end{array}$$(11)

Figure 1 shows F_all as a function of the radius derived from Eq. (11) for three different values of the density ρ₀ of Eq. (8), namely ρ₀ = 5,15,30(× 10⁷) M_⊙ kpc⁻³. The maximum values of the spiral F-strength were 5%, 15%, and 30%, respectively. In the intermediate model, the F-strength varies between 5% and 15% in the region 5 kpc ≤r ≤ 15 kpc. We note that the observed amplitudes (in F−strength) of the spirals in grand-design galaxies with respect to the axisymmetric background have typical values between 5% and 10% (Patsis et al. 1991; Grosbøl et al. 2004).

Fig. 1. Total force perturbation Fall for different amplitudes of the spiral potential as a function of the radius.

Figure 2a shows the surface density profile corresponding to the sum of the disk and bulge component, while Fig. 2b shows the rotation curve produced by the potential V_ax by the equation: $U_{\mathrm{rot}}(r)=\sqrt{r\frac{\partial V_{\mathrm{ax}}(r)}{\partial r}}$.

Fig. 2. (a) Surface density profile of the axisymmetric part of the model. The green curve represents the bulge, the blue curve represents the disk, and the black curve is the total profile. (b) Rotation velocity curve produced by the axisymmetric potential Vax.

Figure 3 shows an isodensity color map of the projected surface density $\sigma (x,\mathit{y})={\int }_{-\infty }^{\infty }\rho (x,\mathit{y},z)\mathrm{d}z$ in the disk plane, where the density ρ is computed from Poisson’s equation ∇²V = 4πGρ for the 3D potential model V(r, ϕ, z) = V_d(r)+V_b(r)+V_sp(r, ϕ, z). It is easy to check that the function ρ(r) is positive everywhere with this choice of potential model (see also Pettitt et al. 2014 and references within and Efthymiopoulos et al. 2020). The maxima of the density are plotted in red. Superimposed are the spiral arms (black curves) derived from the minima of the spiral potential of Eq. (8), which almost completely coincide with the maxima of the projected density of the galactic model. We note that the phase differences between the minima of the potential and the maxima of the spiral density have been suggested to be important close to and outside the corotation region (for a review, see Zhang 2018). However, the spirals considered in our work end at the 4:1 resonance (see below), which is reached well before the corotation region. In the following plots, we therefore use the minima of the potential as indicating the shapes of the imposed spirals for simplicity.

Fig. 3. Isodensities of the projected surface density of the 3D galactic model (see text). Superimposed are the spiral arms derived from the minima of the spiral potential of Eq. (8).

When a fixed pattern speed Ω_sp of the spirals is assumed, the Hamiltonian of stellar orbits in the disk plane in the rotating frame of reference, in polar coordinates, can be expressed as

$$\begin{array}{c}\hfill H=\frac{p_r^2}{2}+\frac{p_{\phi }^2}{2r^2}-{\mathrm{\Omega }}_{\mathrm{sp}}{p}_{\phi }+V_{\mathrm{ax}}(r)+V_{\mathrm{sp}}(r,\phi ),\end{array}$$(12)

where p_r and p_φ are the values of the radial velocity and angular momentum per unit mass in the rest frame.

The angular velocity Ω(r_c) of a star that moves in a circular orbit of radius r_c under the action of the axisymmetric potential alone is given by the relation

$$\begin{array}{c}\hfill \mathrm{\Omega }(r_c)=\sqrt{\frac{1}{r_c}\frac{\mathrm{d}{V}_{\mathrm{ax}}(r_c)}{\mathrm{d}{r}_c}}.\end{array}$$(13)

The epicyclic frequency κ(r_c) at r = r_c is given by

$$\begin{array}{c}\hfill \kappa (r_c)=\sqrt{\frac{{\mathrm{d}}^2{V}_{\mathrm{ax}}(r_c)}{\mathrm{d}{r}_c^2}+\frac{3}{r_c}\frac{\mathrm{d}{V}_{\mathrm{ax}}(r_c)}{\mathrm{d}{r}_c}}.\end{array}$$(14)

Figure 4 shows the function Ω(r), as well as the combinations Ω − κ/2, Ω + κ/2 and Ω − κ/4. The section of Ω(r) with the pattern speed Ω_sp defines the radius of the corotation, while the section of the frequency Ω − κ/2 (or Ω + κ/2) with the pattern speed Ω_sp defines the radius of the ILR (or the outer Lindblad resonance, OLR). Depending on the model, we may have one or two ILRs. Our model (Fig. 4) has two ILRs, the first and the second ILR, as long as Ω_sp < 20 km s⁻¹ kpc⁻¹. Finally, the section of the frequency Ω − κ/4 with the pattern speed Ω_sp defines the radius of the 4:1 resonance. We studied the orbits for three different values of the pattern speed of the spiral potential Ω_sp = 10, 15, and 20 km s⁻¹ kpc⁻¹ (see Sect. 4.2.2).

Fig. 4. Form of the function Ω(r) and of the resonant combinations of Ω(r) and κ(r), Ω − κ/2, Ω + κ/2, and Ω − κ/4. The selected pattern speed Ωsp determines the radii of the ILR, the 4:1 resonance, and the corotation.

3. Extent of the spiral density waves

It is well known that spiral density waves cannot extend throughout the entire galactic disk because natural inner and outer barriers limit the extension of these waves. The spiral density waves are located between the ILR and OLR, but do not reach them (see Dobbs & Baba 2014 for a review). The inner natural barrier approximately coincides with the radius of the 2:1 (or ILR) resonance. The density waves are reflected in this central region before they reach the ILR, and then they are amplified by swing amplification (Goldreich & Tremaine 1978). Furthermore, the density waves can be absorbed at the ILR due to Landau damping (Lynden-Bell 1972). This absorption of the stellar density waves at the ILR can be avoided, however, if the Toomre Q parameter (Toomre 1964) increases significantly (forming a so-called Q-barrier), reflecting the density wave outside the ILR. The Toomre Q parameter for a stellar disk is given by the relation

$$\begin{array}{c}\hfill Q=\frac{\kappa {\sigma }_{\mathrm{R}}}{3.36G{\mathrm{\Sigma }}_0},\end{array}$$(15)

where κ is the epicyclic frequency, σ_R is the velocity dispersion, and Σ₀ is the surface density. An increase in Q parameter signifies high values of the velocity dispersion. This is the case inside the ILR also because of the central spheroidal (bulge). With these assumptions, approximate ‘standing-wave’ patterns can exist between a reflecting radius in the inner part of the galaxy and the corotation radius (Bertin et al. 1989).

However, Contopoulos & Grosbøl (1986, 1988) have shown that when the amplitude of the spiral arms is strong enough, the outer natural barrier coincides with the appearance of the 4:1 (or ultraharmonic) resonance, which is located inside the corotation. The elliptical periodic orbits become rectangular there and can no longer support the spiral density wave beyond that resonance. This result was further confirmed in Patsis et al. (1991, 1994, 1997), Lepine et al. (2011), and Junqueira et al. (2013). Furthermore, Chaves-Velasquez et al. (2019) showed that the spiral arms of a 3D potential model are supported by orbits associated with a stable 2D elliptical periodic orbit as well as its vertical bifurcations. The thickness of the spirals supported by such orbits is compatible with the thickness of the Milky Way disk. However, there are cases where weaker spirals extend beyond the 4:1 resonance (see, e.g., Grosbøl & Patsis 1998).

As an overall conclusion, the spiral density waves supported by precessing ellipses should extend in a region starting from slightly outside the ILR and up to the 4:1 resonance. The radii of these resonances are defined by the specific pattern speed. This region is limited in Fig. 4 between the functions Ω − κ/2 and Ω − κ/4, which is different for different pattern speeds (red line).

4. Phase-space structure

4.1. Periodic orbits

We now describe the main body of our analysis, which is the study of the phase-space structure and of the orbits supporting the spirals in the region between ILR and corotation in the model of Sect. 2. We chose various values for the spiral amplitude (parameter ρ₀ in Eq. (8)), pattern speed Ω_sp, and pitch angle a. Our study focuses on the form and stability of periodic orbits that support the spiral arms as well as the shape of the phase space around these orbits.

Families of the stable periodic orbits that have shapes of precessing ellipses can be found by the Hamiltonian (12), and they correspond to the continuation of the circular orbits of the axisymmetric part of the potential in the region of the 2:1 resonance. The study of such orbits is greatly facilitated using the action-angle variables of epicyclic theory. These are the pair (φ, p_φ) of Eq. (12), as well as the radial angle and action variables (ϑ_r, J_r), defined by

$$\begin{array}{c}\hfill (r-r_c)=\sqrt{\frac{2J_r}{\kappa (r_c)}}sin({\vartheta }_r),p_r=\sqrt{2J_r\kappa (r_c)}cos({\vartheta }_r).\end{array}$$(16)

In Eq. (16), r_c represents a radius of a circle in the disk around which we wish to study the form of the phase portrait. The corresponding Jacobi energy E_j is the energy of the circular orbit of radius r_c under the axisymmetric potential V_ax, H_j(r_c) = H_ax(r_c)−Ω_spp_ϕ(r_c), where $H_{\text{ax}}=p_r^2/2+p_{\varphi }^2/(2r^2)+V_{\text{ax}}$ and $p_{\varphi }(r_c)=r_c^2\mathrm{\Omega }(r_c)$.

Consider, now, the slow angle ψ = ϑ_r − 2φ, as well as the Poincaré canonical variables:

$$\begin{array}{c}\hfill \xi =\sqrt{\frac{2J_r}{\kappa (r_c)}}sin(\psi ),P_{\xi }=\sqrt{2J_r\kappa (r_c)}cos(\psi )\end{array}$$(17)

Using Eqs. (16) and (17), the functions ξ = ξ(r, p_r, φ) and P_ξ = P_ξ(r, p_r, φ) are

$$\begin{array}{c}\hfill \xi =(r-r_c)cos(2\phi )-\frac{p_r}{\kappa (r_c)}sin(2\phi )\end{array}$$(18)

$$\begin{array}{c}\hfill P_{\xi }=p_{r}cos(2\phi )+(r-r_c)\kappa (r_c)sin(2\phi ).\end{array}$$(19)

Equations (18) and (19) can now be used in order to construct a Poincaré surface of section (ξ, P_ξ) for a fixed value of φ and a fixed Jacobi constant E_j. We chose the value φ = π/2, and we find a sequence of surfaces of sections in our model defined by the procedure described below.

For a specific value of the pattern speed Ω_sp, we specified a certain value for the radius of the circular orbit r_c and calculated the corresponding angular momentum $p_{{\varphi }_c}=r_c^2\mathrm{\Omega }(r_c)$ and the corresponding Jacobi integral E_j(r_c) = H_ax(r_c)−Ω_spp_φc. For the value of the total Hamiltonian H = E_j(r_c), we then defined various initial values ξ₀ and P_ξ0, taking as initial value of φ the value φ₀ = π/2. Then we calculated the initial values of the coordinates r, p_r, and p_φ using Eqs. (18) and (19), that is, r₀ = r_c − ξ₀ (from Eq. (18)), p_r0 = −P_ξ0 (from Eq. (19)) and $p_{{\phi }_0}=r_0^2\mathrm{\Omega }(r_0)$. Finally, we used Hamilton’s equations

$$\begin{array}{c}\hfill \frac{\mathrm{d}r}{\mathrm{d}t}=\frac{\vartheta H}{\vartheta p_r},\frac{\mathrm{d}\phi }{\mathrm{d}t}=\frac{\vartheta H}{\vartheta p_{\phi }},\frac{\mathrm{d}{p}_r}{\mathrm{d}t}=-\frac{\vartheta H}{\vartheta r},\frac{\mathrm{d}{p}_{\phi }}{\mathrm{d}t}=-\frac{\vartheta H}{\vartheta \phi },\end{array}$$(20)

with the Hamiltonian (12) in order to integrate orbits with initial conditions (φ₀,r₀,p_r0, and p_φ0) and find the consecutive iterates (ξ, P_ξ) on the Poincaré section φ = π/2.

Figure 5 shows the phase portraits (surfaces of section (ξ, P_ξ)) for the parameter ρ₀ = 5 × 10⁷ M_⊙ kpc⁻³ in Eq. (8) and for the pattern speed Ω_sp = 15 km s⁻¹ kpc⁻¹. This value of the pattern speed places the second ILR not far from the inner break of the surface density profile of the axisymmetric part of the potential due to the presence of the bulge (see Fig. 2a). Figure 5 shows the phase portraits for 12 different values of the radius r_c, namely r_c = 1, 2, ..., 12 kpc, spanning a region from the center of the galaxy and up to a radius just outside the 4:1 resonance.

Fig. 5. Phase-space portraits (ξ, Pξ) for the model of Eq. (8) with pattern speed Ωsp = 15 km s−1 kpc−1 and density of the spiral potential ρ0 = 5 × 107 M⊙ kpc−3 for 12 different values of the radius rc namely rc = 1, 2, ..., 12 kpc. Precessing ellipses responsible for the spiral density waves are periodic orbits of the x1 family that correspond to radii between approximately 5 kpc (second ILR) and 11 kpc (4:1 resonance).

Regarding the main families of periodic orbits in the phase portraits of Fig. 5, we refer to the nomenclature of Contopoulos (1975), who called x₁, x₂ the families of stable and x₃ the family of unstable periodic orbits. In our model, the precessing ellipses supporting the spirals are the stable periodic orbits of the x₁ family. Four different regions can be distinguished according to the number and stability of periodic orbits: (a) Inside the first ILR, only the x₁ family (stable) exists, (b) between the first and the second ILR, three families exist, namely x₁, x₂ (stable), and x₃ (unstable). (c) Between the second ILR and the 4:1 resonance, only the x₁ family exists (stable or unstable family), and (d) outside the 4:1 resonance, the x₁ family still exists, but it no longer supports the spiral arms (see the Figs. 7, 10, 13, and 16). For Ω_sp = 15 km s⁻¹ kpc⁻¹, the second ILR and the 4:1 resonance are approximately at r = 5 kpc and 11 kpc, respectively (Fig. 4). The spirals should thus extend roughly between these two radii.

In the phase portrait of Fig. 5, the x₁ family is found at the center of the islands of stability, marked with black points, while the x₂ family is at the center of the islands, marked with red points. The unstable periodic orbit x₃ is plotted with a blue dot. The orbits x₃ exist only in the region between the first and second ILR, and they do not support the imposed spirals (see below). On the other hand, as shown in Fig. 5, the x₁ family remains stable at all radii up to the 4:1 resonance. Beyond this resonance, orbits of greater multiplicity bifurcate from the x₁ family and substantially affect the structure of the phase space, as in the last panel of Fig. 5.

Figure 6 shows the normalized characteristic curves $S(r_c)=\sqrt{{\xi }^2+P_{\xi }/{\kappa }_c^2}$ (where the epicyclic frequency κ_c(r_c) is given by Eq. (14)) as a function of the radius r_c for the same parameters as in Fig. 5. The family x₁ is shown in black, x₂ in red, and the (unstable) family x₃ in blue. We note that the x₂ and x₃ families exist between the first and second ILR. Both families are created by a tangent bifurcation close to the first ILR, and then rejoin and disappear by inverse bifurcation close to the second ILR. S(r_c) in each case represents the amplitude of the epicyclic oscillation of the corresponding (elliptical) periodic orbit.

Fig. 6. Characteristic curves $S(r_c)=\sqrt{{\xi }^2+p_{\xi }^2/{\kappa }_c^2}$ of the periodic families x1 (black), x2 (red), and x3 (blue) for the model of Eq. (8) with Ωsp = 15 km s−1 kpc−1 and density of the spiral potential ρ0 = 5 × 107 M⊙ kpc−3. The dashed black vertical lines correspond to the first and second ILR, and the dashed gray vertical line corresponds to the 4:1 resonance.

A key remark in Fig. 6 is that in the range of the radii where S(r_c) decreases, between the second ILR (at r_c ≈ 5 kpc) and r_c ≈ 8.5 kpc, the ellipses become more circular, and so they do not intersect with each other. The value of S(r_c) reaches a minimum value of about 8.5 kpc and then increases, first smoothly, and then abruptly. After this latter radius, the orbits become rectangular as they approach the 4:1 resonance. Hence, the x₁ orbits start intersecting and no longer support the imposed spirals. As a conclusion, the end of the response spiral density wave should be placed somewhere between the radius r_c corresponding to the minimum of the curve S(r_c) and the 4:1 resonance.

As shown in Fig. 7, the orbits of the x₁ family form precessing ellipses supporting the imposed spirals in the whole region between the second ILR (≈5 kpc) and a short distance inside the 4:1 resonance. Figure 7a shows that the x₁ family forms a dense and well-defined spiral density wave. In Fig. 7b, the x₁ family is plotted both inside and outside the second ILR. Between the first ILR (≈1.5 kpc) and the second ILR (≈5 kpc), the orbits form a more fuzzy spiral density wave, while inside the first ILR, the orbits are more circular. In fact, inside the second ILR, the spiral density wave is greatly weakened because the amplitude of the spiral perturbation is close to zero (see Fig. 1). The imposed spiral arms (superposed in the figure in red) are derived from the minima of the spiral potential of Eq. (8). The coincidence between imposed spirals and those formed by the elliptical orbits is nearly complete.

Fig. 7. (a) Spiral density waves formed by the precessing ellipses of the elliptical closed orbits of the x1 family, from the model of Eq. (8) for a pattern speed Ωsp = 15 km s−1 kpc−1 and density of the spiral potential ρ0 = 5 × 107 M⊙ kpc−3, between the second ILR and the 4:1 resonance. (b) Same as in (a), but with the ellipses of the x1 family inside the second ILR. Superposed are circles corresponding to the second ILR (green circle), 4:1 resonance (blue circle), and corotation (black circle). The imposed spiral arms (in red) are derived from the minima of the spiral potential of Eq. (8). The coincidence of the imposed spirals with the spiral density wave created by the precessing ellipses is nearly complete.

4.2. Parametric study

The agreement between the imposed spirals (minima of the potential of Eq. (8)) and the response spirals (formed by the elliptical periodic orbits), as well as the regularity of the phase space structure observed in the previous example, holds for a particular choice of parameters (ρ₀, Ω_sp, and a). We now study how the above picture is altered by varying independently the amplitude ρ₀, the pattern speed Ω_sp, or the pitch angle a in the imposed spirals. The results of this parametric study are summarized in Sect. 4.2.3.

4.2.1. Role of the amplitude of the perturbation

In order to investigate the role of the spiral amplitude ρ₀ (Eq. (8)), we repeated the study of the previous subsection for increasing values of ρ₀. Figure 8 shows the same phase-space portraits (ξ, P_ξ) as in Fig. 5, but with ρ₀ three times larger (ρ₀ = 15 × 10⁷ M_⊙ kpc⁻³). When the amplitude is increased, the main observation is that while the x₁ family remains stable for most values of r_c, chaos is introduced in the phase space for radii r_c > 5 kpc covering a great part of the phase space around the island of stability corresponding to the x₁ periodic orbit. The x₁ family itself becomes unstable within a small interval of r_c values (see Fig. 9). Figure 9 shows the normalized characteristic curves $S(r_c)=\sqrt{{\xi }^2+P_{\xi }/{\kappa }_c^2}$ in the model with ρ₀ = 15 × 10⁷ M_⊙ kpc⁻³. As in Fig. 6, here the x₂ and x₃ periodic orbits are also created together at a tangent bifurcation close to the first ILR, and then they join and disappear close to the second ILR. The value of S(r_c) reaches a minimum value around 8.5 kpc (same as in Fig. 6), and then it increases abruptly. Therefore the end of the response spiral density wave here again should be placed somewhere between the radius r_c corresponding to the minimum S(r_c) and the radius corresponding to the 4:1 resonance.

Fig. 8. Same as Fig. 5, but for ρ0 = 15 × 107 M⊙ kpc−3.

Fig. 9. Same as in Fig. 6, but for the density parameter of the spiral potential ρ0 = 15 × 107 M⊙ kpc−3. The dashed part of the black curve denotes the unstable periodic orbit x1.

Figure 10 shows the spiral density waves generated by the x₁ family in this model, extending for radii between the second ILR (r_c ≈ 5 kpc) and a short distance inside the 4:1 resonance. The response spirals now appear more concentrated than the spirals in Fig. 7, around the locus of the maximum of the density wave. This is because the forced ellipticity of the x₁ orbits increases with ρ₀ (see Efthymiopoulos 2010 for a review). In Fig. 10b, the x₁ family is plotted both inside and outside the second ILR. Here again, the orbits between the first and second ILR (from ≈1.5 kpc to ≈5 kpc) form a more fuzzy spiral density wave, while inside the first ILR, the orbits are more circular. On the other hand, we again observe a nearly complete coincidence between imposed and response spirals beyond the second ILR. Figure 11 shows the phase portraits (ξ, P_ξ) for an even greater (by a factor 6) value of ρ₀ = 30 × 10⁷ M_⊙ kpc⁻³. In comparison with Fig. 8, we observe that here chaos is introduced at approximately the same values of r_c (i.e., for r_c > 5 kpc) as in the previous example. However, the main qualitative difference between these two cases is that in the latter case, the x₁ family becomes unstable by a sequence of period-doubling bifurcations starting at r_c ≈ 6.8 kpc, and there are no ordered orbits around it from there on. As a consequence, no ordered matter is collected around it that might support a realistic spiral density wave.

Fig. 10. Same as in Fig. 7, but for ρ0 = 15 × 107 M⊙ kpc−3.

Fig. 11. Same as Fig. 5, but for ρ0 = 30 × 107 M⊙ kpc−3.

Figure 12 shows the normalized characteristic curves $S(r_c)=\sqrt{{\xi }^2+P_{\xi }/{\kappa }_c^2}$ of this model. The main difference with respect to the previous cases is that the curve S(r_c) forms a nearly constant plateau from the second ILR to the point r_c ≈ 7 kpc and then marks an abrupt fall to a minimum at r_c ≈ 8 kpc. As a consequence, the ellipticity of the x₁ orbits remains nearly constant in the region between the second ILR and the radius r_c ≈ 7 kpc. Moreover, when the elliptical orbits of the unstable x₁ periodic orbit are plotted for this model (Fig. 13a), the ellipses intersect each other in the whole range of radii and therefore only define fuzzy spiral density waves (compare with Fig. 10). Some secondary spiral arms also appear. Overall, this is an unrealistic spiral density wave that cannot be observed in real galaxies, as there exists no ordered matter around the x₁ family. Fig. 13b shows the same information as Figs. 7b and 10b for ρ₀ = 30 × 10⁷ M_⊙ kpc⁻³.

Fig. 12. Same as in Fig. 6, but for ρ0 = 30 × 107 M⊙ kpc−3. The dashed part of the black curve denotes that the periodic orbit x1 is unstable.

Fig. 13. Same as in Fig. 7, but for ρ0 = 30 × 107 M⊙ kpc−3.

Comparing all three models, we find that the precessing ellipses of the x₁ family can support the spirals for amplitudes (F-strengths) not exceeding the level 15 − 20%. Beyond this value, the x₁ family becomes unstable at a rather small distance beyond the second ILR (Δr_c ≈ 2 kpc), while chaos dominates the phase space. This is in accordance with estimates of the amplitude of the spiral perturbation of the spiral arms in real grand-design galaxies, which give a relatively low upper limit (≈10 − 15%) in the forces (Grosbøl & Patsis 1998; Grosbøl et al. 2004). Moreover, long-term evolution of self-gravitating models shows that spiral density waves do not remain viable over many revolutions if the spiral forcing is higher than 5% (Chakrabarti et al. 2003).

For the role of other families of periodic orbits, Fig. 14 shows the precessing ellipses of the x₂ family for the model of Eq. (8) for a pattern speed Ω_sp = 15 km s⁻¹ kpc⁻¹ and a density of the spiral potential ρ₀ = 5 × 10⁷ M_⊙ kpc⁻³ in Fig. 14a, ρ₀ = 15 × 10⁷ M_⊙ kpc⁻³ in Fig. 14b, and ρ₀ = 30 × 10⁷ M_⊙ kpc⁻³ in Fig. 14c. In all three cases, the response spirals are orthogonal to the response spirals of the x₁ family (see Contopoulos 1975). In Figs. 14b and c, they form weak spiral arms in the region between the first and second ILR, where the amplitude of the spiral potential is close to zero (see Fig. 1).

Fig. 14. x2 family of orbits for the model of Eq. (8) for (a) a pattern speed Ωsp = 15 km s−1 kpc−1 and density of the spiral potential ρ0 = 5 × 107 M⊙ kpc−3 in (a), ρ0 = 15 × 107 M⊙ kpc−3 in (b) and ρ0 = 30 × 107 M⊙ kpc−3 in (c). In all three cases, these elliptical orbits do not support the spiral wave derived from the galactic potential, and their main axes are perpendicular to the main axes of the x1 family.

In Fig. 15, some precessing ellipses of the x₁ family (black orbits) are plotted together with those of the x₂ family (red orbits) for ρ₀ = 30 × 10⁷ M_⊙ kpc⁻³ and Ω_sp = 15 km s⁻¹ kpc⁻¹, in the range of radii from r_c = 2 kpc to r_c = 5 kpc. The main axes of the ellipses of the x₂ family are perpendicular to the main axes of the x₁ family, and they exist for radii for which the amplitude of the spiral perturbation is close to zero (see Fig. 1). Therefore they do not support the spiral density wave. It can easily be verified that the same is true for the unstable x₃ family (not shown in the figures).

Fig. 15. Some precessing ellipses of the x1 family (black orbits) together with the precessing ellipses of the x2 family (red orbits) that correspond to the same radii rc for density of the spiral potential ρ0 = 30 × 107 M⊙ kpc−3 and for a pattern speed Ωsp = 15 km s−1 kpc−1. The main axes of the x2 ellipses are perpendicular to the main axes of the x1 family and do not support the spiral density wave.

4.2.2. Role of the pattern speed

In order to examine the dependence of the response spiral density wave on the pattern speed Ω_sp, we fixed the amplitude of the spiral perturbation to ρ₀ = 5 × 10⁷ M_⊙ kpc⁻³ and changed the pattern speed of the spiral potential, comparing the cases Ω_sp = 10, 15 and 20 km s⁻¹ kpc⁻¹. Figure 16 shows the spiral density waves formed by the precessing ellipses of the x₁ family (black orbits) for the pattern speed Ω_sp = 20 km s⁻¹ kpc⁻¹ (in Fig. 16a) and Ω_sp = 10 km s⁻¹ kpc⁻¹ (in Fig. 16b). Superposed are the circles corresponding to the ILR radius (green circle), to the 4:1 resonance (blue circle), and to corotation (black circle). By comparing Figs. 16a and b and Figs. 7a and b, which all have the same amplitude ρ₀ = 5 × 10⁷ M_⊙ kpc⁻³, but different values of the pattern speed we derive the following conclusions. When the pattern speed decreases, (a) all the resonances are shifted outward (see Fig. 4), and therefore the spiral density waves reach larger radii. However, the ellipses become rounder when they approach the 4:1 resonance, and therefore the spiral density wave becomes less conspicuous at larger radii. (b) The region inside the first ILR becomes smaller, and the region between the first and second ILR increases. The elliptical orbits of the x₁ family in this latter region become much more elongated and intersect each other. (c) The width of the spiral arms grows with radius (see Fig. 16b). Savchenko et al. (2020) claimed that in the 86% of a sample of 155 face-on grand-design spiral galaxies, they observed that the arm width increases with radius (see Mosenkov et al. 2020, for an alternative interpretation of this phenomenon based on the mechanism of swing amplification).

Fig. 16. Spiral density waves formed by the elliptical orbits of the x1 family, from the model of Eq. (8) for a pattern speed Ωsp = 20 km s−1 kpc−1 in (a) and Ωsp = 10 km s−1 kpc−1 in (b) and a density of the spiral potential ρ0 = 5 × 107 M⊙ kpc−3. Superposed are circles corresponding to ILR (green circle), 4:1 resonance (blue circle), and corotation (black circle).

4.2.3. Role of the pitch angle

In this subsection, we chose the model ρ₀ = 15 × 10⁷ M_⊙ kpc⁻³, Ω_sp = 15 km s⁻¹ kpc⁻¹, but varied the pitch angle from the low value a = −5° to the high value a = −25° instead of the intermediate value a = −13° that was used in all previous numerical experiments.

Figure 17 shows the phase portraits (ξ, P_ξ) for a = −25°. By comparing the two figures (8 and 17, which only differ in the value of the pitch angle), we conclude that for an increasing pitch angle (more open spiral arms), more order is introduced in the phase space and the chaotic areas shrink. Therefore more matter exists in ordered motion around the x₁ stable periodic orbit that can support the spiral density wave better. The corresponding response spirals (Fig. 19) are also more open and extend throughout the whole region between the second ILR and the 4:1 resonance.

Fig. 17. Same as in Fig. 8, but for a pitch angle a = −25°.

Figure 18 instead shows the phase portraits (ξ, P_ξ) when a = −5°. The key observation is that the x₁ family now becomes unstable slightly outside the second ILR (at r_c ≈ 6.65 kpc) and up to the 4:1 resonance, and the region around it is fully chaotic. This case presents an entirely different behavior than all previous cases. We studied the bifurcations of the x₁ and x₂ families in this case in the Appendix A in detail. The corresponding spiral density wave formed by the x₁ family and its continuation, the b family (of Appendix A), is plotted in Fig. 20a. The ellipses of the x₁ family exist up to the second ILR (at r_c ≈ 5 kpc) and intersect, forming a rather fuzzy density wave, while the ellipses from the b family (which is the continuation of the x₁ family) create a more clearly defined density wave that reaches from the second ILR up to the 4:1 resonance (at r_c ≈ 11 kpc). Superposed, in red, are the theoretical minima of the imposed spirals. We observe that they coincide in general with the spiral density wave formed by the precessing ellipses of the x₁ family, although a secondary spiral wave appears from a certain radius and beyond. This is an unrealistic spiral density wave as it is formed by unstable periodic orbits of the x₁ family that have no ordered matter around them, but only chaotic orbits all the way between the second ILR and the 4:1 resonance. Figure 20b instead shows the spiral density wave formed by the x₂ family, which in this model extends far beyond the second ILR (see the Appendix). As commented above, the main axes of the corresponding ellipses are always perpendicular to those of the x₁ family ellipses (see Fig. 15), and thus cannot support the imposed spiral arms.

Fig. 18. Same as in Fig. 8, but for a pitch angle a = −5°.

Fig. 19. Spiral density wave formed by the precessing ellipses of the x1 family for ρ0 = 15 × 107 M⊙ kpc−3, a pattern speed Ωsp = 15 km s−1 kpc−1, and a pitch angle a = −25°.

Fig. 20. Spiral density waves formed by the precessing ellipses of the (a) x1 family and (b) x2 family for ρ0 = 15 × 107 M⊙ kpc−3, a pattern speed Ωsp = 15 km s−1 kpc−1, and a pitch angle a = −5°. Superposed (in red) is the theoretical spiral derived from the minima of the spiral potential of Eq. (8).

By comparing Figs. 8, 17, and 18, we conclude that for an increasing pitch angle (more open spiral arms), more order is introduced in the phase space. For very small pitch angles, chaos dominates in the phase space, and we always obtain that the x₁ family of orbits becomes unstable almost immediately after the second ILR. The x₂ family is also unstable in the region outside the second ILR.

In conclusion, there is a lower limit of the value of the pitch angle a below which the central periodic orbit of the x₁ family becomes unstable in the region between the second ILR and the 4:1 resonance. This limit of the pitch angle is defined by the other two free parameters of the model, that is, by the amplitude of the spiral perturbation and the pattern speed of the spiral arm, as indicated in Table 1.

Table 1 can be used to estimate the permissible area of pitch angles that for ρ₀ and Ω_sp as in the first two columns should be larger than the value reported in the third column as a function of the amplitude of the spirals and the pattern speed. From these data, we conclude that using the stability of the x₁ family as a criterion, a correlation between the pitch angle and the amplitude of the spiral perturbation suggests that spirals formed by precessing ellipses should be stronger in amplitude when they are more open (larger a). Moreover, a correlation between the pitch angle and the pattern speed is indicated: the higher the value of the pattern speed, the more tightly wound the spirals should be to maintain stability of the x₁ family.

4.3. Comparison with observations

In its classical version, the Hubble sequence for spiral galaxies implies that the bulge size and the spiral arm winding should be highly correlated. According to this classification, the “Sa” galaxies have tightly wound arms and fat nuclear bulges, “Sb” galaxies have moderately wound arms and moderate nuclear bulges, and “Sc” galaxies have loosely wound arms and small nuclear bulges. However, modern classification schemes for spiral galaxies imply a considerable departure from the classic Hubble sequence as regards the correlation between spiral arm winding type and bulge size. Early studies suggested that spiral arms become tighter with increasing mass of the bulge (Morgan 1958, 1959; Kennicutt 1981; Bertin et al. 1989). Furthermore, Seigar et al. (2005, 2006) reported a tight connection between pitch angle and morphology of the galactic rotation curve, quantified by the shear rate, with open arms associated with rising rotation curves and tightly wound arms connected to flat and falling rotation curves. On the other hand, Hart et al. (2017) recently analyzed a large sample of galaxies selected from the Sloan Digital Sky Survey (SDSS; York et al. 2000) and found very weak correlations between pitch angle and galaxy mass, and a surprising trend that the pitch angle increases with increasing bulge-to-total mass ratio. Yu & Ho (2019) found that the pitch angle decreases (arms are more tightly wound) in galaxies of earlier Hubble type, more prominent bulges, higher concentration, and higher total galaxy stellar mass. However, there is a significant scatter in their measures. Finally, Masters et al. (2019) and Díaz-García et al. (2019) found little or no correlation between spiral arm winding tightness and bulge size.

There is no concluding evidence of a correlation between the size of the bulge and the pitch angle of the grand design galaxies. We therefore tested various pitch angles in our galactic model and kept the mass of the bulge constant at a relatively high value. The results of the previous subsection indeed suggest a correlation between pitch angle and amplitude of the spiral perturbation, suggesting that galaxies with stronger perturbations should have more open spirals (greater value of the pitch angle). In comparison with real observations, Grosbøl et al. (2002) found that the distribution of mean amplitudes of two-armed spirals as a function of the pitch angle shows a lack of strong in amplitude and at the same time tightly wound spirals. They have also found that most of the mean amplitudes of spiral arms are below 15% in forces, and there is a correlation between the amplitude of spirals and the pitch angle: weak spirals often have tighter spiral arms. Díaz-García et al. (2019) also found that the mean amplitude of the arms increases with increasing pitch angle.

5. Conclusions

We studied the precessing ellipse model of elliptical orbits that support the spiral density waves in grand-design spiral galaxies. We used a theoretical model consisting of a bulge, a disk, a halo, and a spiral potential. The free parameters of our model were the amplitude of the spiral density perturbation (ρ₀ in Eq. (8)), the pattern speed of the spiral potential (Ω_sp in Eq. (12)), and the pitch angle of the spiral arms of the model (a in Eq. (8)). By testing the effect of the variation in free parameters of the model on the response spirals, formed by elliptical periodic orbits, we extracted the conclusions listed below.

• (1)

  In all models we studied, the x₁ family creates response spirals that are consistent with the imposed spirals. The response spirals extend in the region from the center of the galaxy up to the radius of the 4:1 resonance. The x₂ and x₃ family of orbits exist between the first and second ILR. They are created simultaneously near the first ILR at a tangent bifurcation, and they join and disappear near the second ILR in all the models, except in the case of a very small pitch angle, where they still exist outside the second ILR, but do not contribute to the response spirals. The main axes of the x₂ family of orbits are perpendicular to the main axes of the x₁ family, and the x₃ family is always unstable in the whole range of radii.

• (2)

  By increasing the amplitude of the spiral density perturbation in our model, chaos is introduced gradually and the x₁ family of orbits becomes unstable. An upper limit for the spiral density perturbation is ρ₀ = 30 × 10⁷ M_⊙ kpc⁻³, where the x₁ family is unstable in the whole range between the second ILR and the 4:1 resonance.

• (3)

  When the value of the pattern speed Ω_sp decreases, it causes an outward shift of the resonances and therefore the spiral density waves can reach larger radii. However, the ellipses become rounder when they approach the 4:1 resonance, and therefore the spiral density wave becomes less conspicuous at larger radii. Moreover, the elliptical orbits of the x₁ family become much more elongated and intersect, thus destroying a coherent spiral response.

• (4)

  We conclude for the value of the pitch angle that for an increasing pitch angle (more open spiral arms), more order is introduced in the phase space and the chaotic areas shrink, while for a decreasing pitch angle (more tight spiral arms), chaos dominates in the phase space and the x₁ family of orbits can no longer support the spiral density wave. In the special case of a very small pitch angle, new bifurcations of periodic orbits appear from the main families x₁ and x₂, but the main bifurcations of the x₁ family that cause the spiral density wave become unstable soon after the second ILR, and therefore they cannot support real spirals.

• (5)

  To summarize, the main result of the paper is a quantitative estimation of the range of the free parameters of our galactic model for generating realistic spiral density waves through the precessing ellipses model of elliptical closed orbits (and their surrounding quasi-periodic ones). In particular, a correlation between the pitch angle and the amplitude of the spiral perturbation was shown, where stronger spirals can be statistically less tight than weaker spirals. Moreover, a correlation between the pitch angle and the pattern speed was shown, where spirals that spin faster can be statistically tighter than slower spirals.

Acknowledgments

We acknowledge support by the research committee of the Academy of Athens through the project 200/895.

References

Appendix A: Bifurcations of the orbits in the case of a small pitch angle

Section 4.2.3 describes the case of a very small pitch angle, a = −5°, for the amplitude of the spiral perturbation ρ₀ = 15 × 10⁷ M_⊙ kpc⁻³ and a pattern speed Ω_sp = 15 km s⁻¹ kpc⁻¹ in our galactic model. In this case, the families of periodic orbits emerge from an intricate sequence of bifurcations. In this appendix, we present some further details of these families.

Figure 18 shows the phase-space portraits of these families for integer values of the radius of the circular orbit r_c, namely r_c = 1, 2, 3, ..12 kpc. The various bifurcations appear at values of r_c between those shown in the figure.

Figure A.1 shows the characteristic curves $S(r_c)=\sqrt{{\xi }^2+p_{\xi }^2/{\kappa }_c^2}$ of the various families for this case. Fig. A.1a shows the characteristic curve of the x₁ family and its bifurcations, while Fig. A.1b shows the characteristic curve of the x₂ family and its bifurcations, as well as the x₃ family. Figure A.2 shows the phase-space portraits for radii r_c = 4.4, 4.5, 4.6, 4.7, 4.78, 4.8, 5.0, and 5.5 kpc, which give a more detailed description of Fig. 18 in the range of radii between r_c = 4 and r_c = 5.5, where the bifurcations of the x₁ and x₂ families take place.

Fig. A.1. (a) Characteristic curves $S(r_c)=\sqrt{{\xi }^2+p_{\xi }^2/{\kappa }_c^2}$ of the periodic families x1 (black) and their bifurcations $x_1^{\prime }$ (gray), $x_1^{″}$ (orange), g (green), and b (blue). The dashed part of the blue curve denotes that the b family is unstable. (b) The characteristic curves S(rc) of the periodic families x3 (dashed blue) and x2 (red), and their bifurcations x4 (pink) and x5 (dashed red). The parameters of the model of Eq. (8) are Ωsp = 15 km s−1 kpc−1, ρ0 = 15 × 107 M⊙ kpc−3, and a = −5°.

These figures show that a) for values of r_c close to zero, only the stable x₁ family exists. This family exists all the way up to r_c ≈ 4.84, where it joins the unstable family $x_1^{″}$ and disappears. b) The families x₂ (stable) and x₃ (unstable) are generated at a tangent bifurcation at r_c ≈ 1.2. The family x₂ exists all the way up to r_c ≈ 8.5, where it joins the unstable family x₅ and disappears. The family x₃ exists up to r_c ≈ 4.98, where it joins the unstable family x₄ and disappears. c) The phase-space portrait (ξ, P_ξ) for r_c = 4.4 (Fig. A.2a) shows x₁, x₂, x₃, and some new bifurcated orbits inside the domain covered by red invariant curves, surrounding the x₂ periodic orbit. These orbits are called b (blue stable periodic orbit) and g (green unstable periodic orbit) and have been generated at r_c ≈ 4.23 at a tangent bifurcation (Fig. A.1a). Fig. A.2a also clearly shows that the separatrix from x₃ surrounds x₂ (inner curve of separatrix) and both x₁ and x₂ (outer curve of separatrix). d) From Fig. A.2b (r_c = 4.5), we conclude that for a slightly smaller r_c, the unstable periodic orbit g has moved downward and has reached the separatrix emanating from x₃. e) Fig. A.2c (r_c = 4.6) shows that the separatrix of g surrounds both the x₁ and the b periodic orbits. The separatrix of x₃ no longer surrounds the orbit b. f) For r_c = 4.6 (Fig. A.2c), two new orbits appear in the domain around the x₂ periodic orbit, namely orbits x₄ (stable) and x₅ (unstable) created at a tangent bifurcation at r_c ≈ 4.51 (Fig. A.1b). The separatrix emanating from x₅ surrounds both x₂ and x₄. Moreover, the separatrix of x₃ surrounds all the three orbits, namely x₂, x₄, and x₅. g) Fig. A.2d (r_c = 4.7) shows that two more families of orbits are created inside the region surrounding the x₁ periodic orbit, namely $x_1^{\prime }$ (stable) and $x_1^{″}$ (unstable). These families are generated at a tangent bifurcation at r_c ≈ 4.64 (Fig. A.1a). h) Fig. A.2e (r_c = 4.78) shows that the separatrix emanating from x₃ has almost reached the unstable orbit x₅ and the separatrix emanating from x₅ (that surrounds x₄) has almost reached x₃. For a slightly smaller value of r_c, these parts of the two separatrices should coincide. i) At r_c = 4.8 (Fig. A.2f), the separatrix emanating from x₃ surrounds only x₄. On the other hand, the separatrix emanating from x₅ surrounds on the one side the periodic orbit x₂ and the other side all the periodic orbits $x_1,x_1^{\prime },x_1^{″},b$, and g. j) At r_c ≈ 4.98, the orbits x₃ and x₄ join and disappear, and the same happens for orbits x₁ and $x_1^{″}$ at r_c ≈ 4.84 (Figs. A.1a,b). Thus, for r_c = 5 (Fig. A.2g), only the periodic orbits $x_1^{\prime }$ (stable), b (stable), g (stable), x₂ (stable), and x₅ (unstable) exist. k) At r_c ≈ 5.45, the orbits $x_1^{\prime }$ and g join and vanish (Fig. A.1a). Thus, for r_c = 5.5 (Fig. A.2h), only orbits b (stable), x₂ (stable), and x₅ (unstable) remain. l) At r_c ≈ 8.5, families x₂ and x₅ join and vanish (Fig. A.1b), and finally, beyond r_c ≈ 8.5 only the b periodic orbit remains up to the 4:1 resonance. This is the main family that could support the spiral density wave, but it is unstable in most of the domain between the second ILR and the 4:1 resonance. Although we have constructed a density wave out of the precessing ellipses of the unstable periodic orbit b (see Fig. 20a), the fact that there are no ordered orbits around this periodic orbit, has as a result that such a spiral density wave cannot be observed in real galaxies.

Fig. A.2. Phase-space portraits (ξ, Pξ) for the model of Eq. (8) with pattern speed Ωsp = 15 km s−1 kpc−1, ρ0 = 15 × 107 M⊙ kpc−3, and a = −5°, for values of the radius rc = 4.4, 4.5, 4.6, 4.7, 4.78, 4,8, 5.0, and 5.5 kpc. Precessing ellipses responsible for the spiral density waves are periodic orbits of the b family (bifurcation of the x1 family).

The sequence of bifurcations that we observe in Figs. A.1a and b is much more complicated than in the corresponding Figs. 6, 9, and 12, where only families x₁, x₂ and x₃ exist. In these cases, families x₂ and x₃ start at a tangent bifurcation and disappear by joining again. Family x₁ is stable all the way to the 4:1 resonance and beyond it for a small amplitude of the spiral perturbation (ρ₀ = 5 × 10⁷ M_⊙ kpc⁻³, Fig. 5), but it becomes unstable at a certain radius between the second ILR and the 4:1 resonance by a sequence of period-doubling bifurcations for higher values of the amplitudes, that is, for ρ₀ ⩾ 15 × 10⁷ M_⊙ kpc⁻³ (see Figs. 8 and 11).

In all these cases, the pitch angle is a = −13°. For larger pitch angles (see Fig. 16 for a = −25°), the x₁ family is stable at least up to the 4:1 resonance.

In the present case, where the pitch angle is significantly smaller, that is, a = −5°, family x₁ has a very different evolution. The characteristic curve of x₁ (stable) is continued by the characteristic curve of $x_1^{″}$ (unstable and backward in r_c), then $x_1^{\prime }$ (stable and forward in r_c), g (unstable and backward in r_c), and b (stable and forward in r_c) (Fig. A.1a), which finally becomes unstable by a period-doubling bifurcation slightly after r_c ≈ 6 and produces great chaos. On the other hand, the characteristic curve of x₃ (unstable) starts at tangent bifurcation with x₂ (stable) and continues with the characteristic curve of x₄ (stable, backward in r_c), the x₅ (unstable), and finally, it joins x₂ (stable) at r_c ≈ 8.5 and disappears. Thus the extended characteristic curves of x₃ and x₂ have the same starting and ending points as in all the previous models.

Another strange new feature in this model is the evolution of the separatrices of the new unstable families g and x₅. In both cases, the unstable points of g and x₅ reach the separatrix emanating from the unstable point x₃ at critical values of r_c = r_ccrit. After this value of r_c, the separatrices emanating from g or x₅ change drastically.

All Tables

All Figures

	Fig. 1. Total force perturbation Fall for different amplitudes of the spiral potential as a function of the radius.
In the text

	Fig. 2. (a) Surface density profile of the axisymmetric part of the model. The green curve represents the bulge, the blue curve represents the disk, and the black curve is the total profile. (b) Rotation velocity curve produced by the axisymmetric potential Vax.
In the text

	Fig. 3. Isodensities of the projected surface density of the 3D galactic model (see text). Superimposed are the spiral arms derived from the minima of the spiral potential of Eq. (8).
In the text

	Fig. 4. Form of the function Ω(r) and of the resonant combinations of Ω(r) and κ(r), Ω − κ/2, Ω + κ/2, and Ω − κ/4. The selected pattern speed Ωsp determines the radii of the ILR, the 4:1 resonance, and the corotation.
In the text

Fig. 5. Phase-space portraits (ξ, Pξ) for the model of Eq. (8) with pattern speed Ωsp = 15 km s−1 kpc−1 and density of the spiral potential ρ0 = 5 × 107 M⊙ kpc−3 for 12 different values of the radius rc namely rc = 1, 2, ..., 12 kpc. Precessing ellipses responsible for the spiral density waves are periodic orbits of the x1 family that correspond to radii between approximately 5 kpc (second ILR) and 11 kpc (4:1 resonance).

In the text

	Fig. 6. Characteristic curves $S(r_c)=\sqrt{{\xi }^2+p_{\xi }^2/{\kappa }_c^2}$ of the periodic families x1 (black), x2 (red), and x3 (blue) for the model of Eq. (8) with Ωsp = 15 km s−1 kpc−1 and density of the spiral potential ρ0 = 5 × 107 M⊙ kpc−3. The dashed black vertical lines correspond to the first and second ILR, and the dashed gray vertical line corresponds to the 4:1 resonance.
In the text

Fig. 7. (a) Spiral density waves formed by the precessing ellipses of the elliptical closed orbits of the x1 family, from the model of Eq. (8) for a pattern speed Ωsp = 15 km s−1 kpc−1 and density of the spiral potential ρ0 = 5 × 107 M⊙ kpc−3, between the second ILR and the 4:1 resonance. (b) Same as in (a), but with the ellipses of the x1 family inside the second ILR. Superposed are circles corresponding to the second ILR (green circle), 4:1 resonance (blue circle), and corotation (black circle). The imposed spiral arms (in red) are derived from the minima of the spiral potential of Eq. (8). The coincidence of the imposed spirals with the spiral density wave created by the precessing ellipses is nearly complete.

In the text

	Fig. 8. Same as Fig. 5, but for ρ0 = 15 × 107 M⊙ kpc−3.
In the text

	Fig. 9. Same as in Fig. 6, but for the density parameter of the spiral potential ρ0 = 15 × 107 M⊙ kpc−3. The dashed part of the black curve denotes the unstable periodic orbit x1.
In the text

	Fig. 10. Same as in Fig. 7, but for ρ0 = 15 × 107 M⊙ kpc−3.
In the text

	Fig. 11. Same as Fig. 5, but for ρ0 = 30 × 107 M⊙ kpc−3.
In the text

	Fig. 12. Same as in Fig. 6, but for ρ0 = 30 × 107 M⊙ kpc−3. The dashed part of the black curve denotes that the periodic orbit x1 is unstable.
In the text

	Fig. 13. Same as in Fig. 7, but for ρ0 = 30 × 107 M⊙ kpc−3.
In the text

Fig. 14. x2 family of orbits for the model of Eq. (8) for (a) a pattern speed Ωsp = 15 km s−1 kpc−1 and density of the spiral potential ρ0 = 5 × 107 M⊙ kpc−3 in (a), ρ0 = 15 × 107 M⊙ kpc−3 in (b) and ρ0 = 30 × 107 M⊙ kpc−3 in (c). In all three cases, these elliptical orbits do not support the spiral wave derived from the galactic potential, and their main axes are perpendicular to the main axes of the x1 family.

In the text

Fig. 15. Some precessing ellipses of the x1 family (black orbits) together with the precessing ellipses of the x2 family (red orbits) that correspond to the same radii rc for density of the spiral potential ρ0 = 30 × 107 M⊙ kpc−3 and for a pattern speed Ωsp = 15 km s−1 kpc−1. The main axes of the x2 ellipses are perpendicular to the main axes of the x1 family and do not support the spiral density wave.

In the text

	Fig. 16. Spiral density waves formed by the elliptical orbits of the x1 family, from the model of Eq. (8) for a pattern speed Ωsp = 20 km s−1 kpc−1 in (a) and Ωsp = 10 km s−1 kpc−1 in (b) and a density of the spiral potential ρ0 = 5 × 107 M⊙ kpc−3. Superposed are circles corresponding to ILR (green circle), 4:1 resonance (blue circle), and corotation (black circle).
In the text

	Fig. 17. Same as in Fig. 8, but for a pitch angle a = −25°.
In the text

	Fig. 18. Same as in Fig. 8, but for a pitch angle a = −5°.
In the text

	Fig. 19. Spiral density wave formed by the precessing ellipses of the x1 family for ρ0 = 15 × 107 M⊙ kpc−3, a pattern speed Ωsp = 15 km s−1 kpc−1, and a pitch angle a = −25°.
In the text

	Fig. 20. Spiral density waves formed by the precessing ellipses of the (a) x1 family and (b) x2 family for ρ0 = 15 × 107 M⊙ kpc−3, a pattern speed Ωsp = 15 km s−1 kpc−1, and a pitch angle a = −5°. Superposed (in red) is the theoretical spiral derived from the minima of the spiral potential of Eq. (8).
In the text

Fig. A.1. (a) Characteristic curves $S(r_c)=\sqrt{{\xi }^2+p_{\xi }^2/{\kappa }_c^2}$ of the periodic families x1 (black) and their bifurcations $x_1^{\prime }$ (gray), $x_1^{″}$ (orange), g (green), and b (blue). The dashed part of the blue curve denotes that the b family is unstable. (b) The characteristic curves S(rc) of the periodic families x3 (dashed blue) and x2 (red), and their bifurcations x4 (pink) and x5 (dashed red). The parameters of the model of Eq. (8) are Ωsp = 15 km s−1 kpc−1, ρ0 = 15 × 107 M⊙ kpc−3, and a = −5°.

In the text

	Fig. A.2. Phase-space portraits (ξ, Pξ) for the model of Eq. (8) with pattern speed Ωsp = 15 km s−1 kpc−1, ρ0 = 15 × 107 M⊙ kpc−3, and a = −5°, for values of the radius rc = 4.4, 4.5, 4.6, 4.7, 4.78, 4,8, 5.0, and 5.5 kpc. Precessing ellipses responsible for the spiral density waves are periodic orbits of the b family (bifurcation of the x1 family).
In the text