Open Access
Issue
A&A
Volume 698, June 2025
Article Number L26
Number of page(s) 6
Section Letters to the Editor
DOI https://doi.org/10.1051/0004-6361/202554837
Published online 17 June 2025

© The Authors 2025

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.

1. Introduction

Recently, there has been a flash of interest in unraveling the nature of buckling, the sudden loss of vertical symmetry by a galactic bar (Collier 2020; Sellwood & Gerhard 2020; Li et al. 2023; Łokas 2019, 2025). This violent buckling episode is observed in many numerical models of disk galaxies (e.g., Sotnikova & Rodionov 2003; Martinez-Valpuesta et al. 2006; Smirnov & Sotnikova 2019; Łokas 2019; Collier 2020; Kataria 2024) and often precedes a rapid thickening of the bar and the formation of the so-called boxy- or peanut-shaped (BP) bulge (Combes & Sanders 1981).

In numerical simulations of barred disk galaxies, buckling was first noticed by Friedli & Pfenniger (1990) and then described as a rapid asymmetric protrusion of the bar in numerical simulations of Raha et al. (1991). Raha et al. (1991) linked it to the fire-hose instability of a thin anisotropic layer. At the same time, Combes et al. (1990), Pfenniger & Friedli (1991), Quillen (2002) explained the general thickening of the bar by the lifting of orbits captured by or passing through the emerging vertical resonance. Quillen et al. (2014) constructed first- and second-order Hamiltonian models of vertical resonance and associated bar asymmetry during its rapid thickening with a rapid distortion of the phase space θres − Jz, in which fixed points (corresponding to 3D periodic orbits) for the same value of the Hamiltonian H have different values of Jz for θres = 0 and π (a first order resonance).

Despite the existing investigations concerning the transformation of individual orbits during buckling (e.g., Martinez-Valpuesta et al. 2006; Łokas 2019; Valencia-Enríquez et al. 2023), there is still no consensus on the nature of buckling in full-fledged N-body simulations. Sellwood & Gerhard (2020) insist on fire-hose instability. Li et al. (2023), Łokas (2025) argue for the resonant nature of vertical heating. They focus on the fact that, shortly before buckling, there is a distortion of the Laplace plane in the same direction along the z-axis at both ends of the bar. As a result, a periodic force with a frequency of 2(Ω − Ωp) (where Ω is the angular velocity and Ωp is the bar pattern speed) begins to act on the almost flat orbits that oscillate with a high vertical frequency. As a consequence, the orbit undergoes forced oscillations with a vertical frequency of ωz ≈ 2(Ω − Ωp). In other words, it enters a vertical resonance. Since the perturbing force pushes the orbits at the ends of the bar predominantly in one direction, a vertical asymmetry arises.

In this Letter we present a different picture of what occurs during buckling. We associate the symmetric exit of the bar orbits from the midplane with perturbations of ∼cos(2θres) (the second-order resonance, Quillen et al. 2014), which arise as a consequence of bar formation and its quiet evolution, and create two fixed points (at θres = 0 and π) in the phase space. However, the culprit of the buckling is asymmetric perturbations, ∼cos θres, that distort the phase space. If they are strong enough, they allow only one fixed point in the phase space1 (the first-order resonance, Quillen et al. 2014). Using the framework of angle-action analysis developed for the self-consistent N-body models in our previous papers Zozulia et al. (2024a,b), we show for the first time how individual orbits interact with asymmetric phase-space perturbations in the system, and how these interactions result in changes to the overall demographics of orbits. Thus, we present here the first detailed description of buckling in terms of individual orbit types differentiated based on the behavior of their resonant angle.

2. Numerical model

We used the same numerical model as in Zozulia et al. (2024a,b) and refer to the cited papers for more details. The model initially consisted of two components, a pure stellar disk (exponential) and an isotropic spherical dark halo that has a Navarro-Frenk-White-like density profile (Navarro et al. 1996). The model parameters were presented in the natural system of units, where the gravitational constant G = 1, the mass of the disk Md = 1, and the initial exponential scale of the disk Rd = 1. One time unit would correspond to 13.8 Myr if the total mass of the disk were Md = 5 × 1010M and the scale length were Rd = 3.5 kpc. In dimensionless units the initial thickness of the disk is zd = 0.05, where zd enters in the sech2″ law, the dark halo mass is Mh(R < 4Rd)≈1.5, and the minimum value of the Toomre parameter is Q(2Rd) = 1.2.

The components evolve in a self-consistent manner until t = 600 (∼8 Gyr). A strong and fast bar appears in the center of the disk at about t = 80. Almost immediately after its formation, the bar begins to thicken. In the time interval from t = 160 to t = 190, the bar buckles rapidly and asymmetrically, reaching maximum asymmetry in the vertical direction at t = 180. By this time, a clearly visible BP bulge develops, which after t = 200 gradually becomes symmetric. Subsequently, the flat bar and the BP bulge continue to increase in size, capturing new disk orbits and lifting them out of the disk midplane (Zozulia et al. 2024a,b).

3. Action-angle variables and orbital types

We associated different types of orbits with different areas in the phase space based on the behavior of the resonant angle2: θres = θz − θR. We calculated the evolution of the non-perturbed, medium-term action-angle variables (JR, Jz, Lz, θR, θz, θϕ) directly in the N-body model, as described in Zozulia et al. (2024a,b). These medium-term variables change smoothly on the timescale of short-term orbital oscillations, allowing orbital libration and circulation to be traced. The actions preserve the mean value of instantaneous actions (calculated in an axisymmetric potential using the AGAMA package, Vasiliev 2019) between apocenters or z-maxima. Angles and frequencies were calculated consistently such that θR changes by 2π between apocenters, θz changes by 2π between z-maxima, and θϕ changes by the azimuthal angle traversed by the particle between the apocenters. We also used secular values of Jz, which are obtained by averaging medium-term values, and thus tracked the orbital evolution on libration timescales. Across this Letter we specify which Jz value is used (medium-term or secular).

For each bar particle, we trace the evolution of the resonant angle θres = θz − θR forward and backward in time, identify four possible modes of its behavior, and thus obtain four types of orbits. These four types of orbits may be described as follows:

(vCIR+) Circulation of the resonant angle with positive angular velocity, θ ˙ z > θ ˙ R $ \dot{\theta}_z > \dot{\theta}_{\mathrm{R}} $ (ωz − κ > 0).

(vCIR−) Circulation of the resonant angle with a negative angular velocity, θ ˙ z < θ ˙ R $ \dot{\theta}_z < \dot{\theta}_{\mathrm{R}} $ (ωz − κ < 0). This orbit circulates, if θres passes through at least 2π in either the positive or negative direction. The notations “up” and “down” of both circulations correspond to the value of the resonant angle in the range −π/2 < θres < π/2 and π/2 < θres < 3π/2, respectively.

(BAN) Orbits in a vertical inner Lindblad resonance (vILR), whose θres undergoes at least one libration (revolution) around a stable fixed point near 0 or π (BAN up and BAN down in our notation). Such orbits are typically banana-shaped.

(vPAS) Orbits in the process of changing the direction of circulation pass close to the vILR, but do not get stuck in it. The passage, in turn, can occur near 0 (vPAS up) or near π (vPAS down).

A detailed description of the calculation of the action-angle variables and the determination of different orbital types (including ILR orbits) based on the changes in resonant angle is given in Appendix A of Zozulia et al. (2024b) and in Appendix A of this Letter.

4. Demographics of orbits

Figure 1 shows how the orbital types in the bar change in relation to θres over the time interval from t = 160 to t = 180 (moment of the highest vertical asymmetry of the bar). The amplitude of the bar during buckling in our model practically does not change (see Fig. 3, left plot in Smirnov & Sotnikova 2018, green line). The fraction of particles in the bar of all disk particles from t = 150 to t = 200 changes slightly, from 32.3% to 33.7%. It then slowly increases to 42.1% by t = 300. The vertical structure of the bar at t = 160 is almost symmetric (Fig. 1, ILR). Most of the bar’s orbits are nearly flat (vCIR+, 78.3%), but among these orbits there is a bias toward θres = 0. This bias is compensated for by the reverse bias of the BAN orbits and the vPAS orbits toward θres = π (“down” types).

thumbnail Fig. 1.

Dynamic decomposition of the bar (ILR) at the beginning of buckling (t = 160, two top row) and at the peak of buckling (t = 180, two bottom row). For each time moment, the five leftmost snapshots on the xy-plane, covering [ − 2.5, 2.5]×[−2.5, 2.5], and on the xz-plane, covering [ − 2.5, 2.5]×[−1.5, 1.5], are shown. From left to right: All bar orbits (ILR); orbits with increasing resonant angle, θres = θz − θR (vCIR+); orbits in the vILR (BAN); orbits passing through it (vPAS); and orbits with decreasing θres (vCIR−). The next five snapshots show the same orbits on the xz-plane, but separated into “up” (near θres = 0) and “down” (near θres = π) types, and into secular Jz < 0.02 (very flat) and Jz > 0.02 (thick) for BAN orbits. For more information on the classification of these types, see Sect. 3 and Appendix A.

The picture changes dramatically at t = 180. By this time, almost two thirds of the vCIR+ orbits are transformed into other types, and although the vCIR+ orbits are still dominated by orbits near θres = 0 (18.8% of the total orbits in the bar versus 7.4% of the vCIR+ orbits near θres = π), a strong vertical asymmetry is visible in the combined snapshot. The asymmetry is due to the large percentage of the BAN down and vPAS down orbits near π (36.3% of the total number of the orbits in the bar). These orbits have significantly greater Jz than the vCIR+ orbits. An additional vertical skew is associated with the appearance of very flat (secular Jz < 0.02) BAN up orbits (θres = 0) that are greatly extended along R. Although the total number of all types of orbits near θres = 0 (48.3%) is almost the same as that near θres = π (51.7%), orbits near θres = π are dominated by those that are in the process of transformation and have high values of Jz. This determines the vertical asymmetry.

Figure 2 shows all of the noted features in changes in the percentages of orbits of different types. The sharp bias of the PAS down orbits with respect to the PAS up ones shortly before buckling and the sharp increase in all the PAS orbits at the moment of buckling are especially visible. Additionally, we note that by t = 300, the symmetry between the orbits near θres = 0 and π is practically restored for all types of orbits, although a slight asymmetry still remains.

thumbnail Fig. 2.

The evolution of the fraction between different types of orbits from t = 100 to t = 300. In all plots, the buckling stage (t = 150 − 200) is highlighted by vertical lines. From left to right: among orbits in ILR (the bar), among ILR orbits captured in vILR (BAN), among ILR orbits passing through vILR (heating), among ILR orbits with increasing resonant angle θz − θR (vCIR+, flat orbits) and among ILR orbits with decreasing resonant angle θz − θR (vCIR−). The different types of orbits in plots, starting with the second, are shown in the same colors used in the first plot. All types correspond to those shown in Fig. 1.

5. Phase-space distortion and portals for transferring orbits

In Zozulia et al. (2024b) we presented typical phase trajectories of individual orbits at late stages of bar evolution, after t = 400. Here, to explain the bias in the number of orbits of different types near θres = 0 and π, we construct phase portraits of orbital ensembles during buckling, which give an idea of the structural evolution and distortions of the phase space.

We developed a method for obtaining sketches of phase portraits on the planes (θres = θz − θR,  Jz) and ( J z cos θ res , J z sin θ res ) $ (\sqrt{J_z}\,\cos{\theta_{\mathrm{res}}},\,\sqrt{J_z}\,\sin{\theta_{\mathrm{res}}}) $ and tracked their evolution directly in the N-body model. To do this, we found orbits that align along the major axis of the bar and have the same fixed value of the Jacobi integral, HJ, and track the evolution of their θres and Jz. This approach enabled the reduction of the 6D phase space of the bar to a 2D one, thus facilitating a more detailed study. Here, we used medium-term Jz and θres so as not to miss the libration of the orbits. More detailed information on the construction of phase portraits can be found in Appendix B.

Figure 3 shows the evolution of the phase portraits described above for different values of HJ (higher modulus values correspond to more central orbits) during and after buckling. The layered structure of the phase space can be observed on the plane (θres, Jz) at each moment of time and for all HJ. We can see that at t = 160 vCIR+ (the bottom orbit layer) and BAN down (the upper layer around θres = π) orbits occupy most of the phase space. The area size of the BAN up orbits is negligibly small. For given values of HJ, there is not even a stable fixed point near θres = 0. The areas occupied by vCIR− orbits (the uppermost layer), flat BAN up orbits (the new lowest layer), and non-flat BAN up (the central layer around θres = 0) increase significantly by t = 180 and 200. All of this fits perfectly with the picture that we saw after performing the dynamic decomposition of the bar (Fig. 1). Also, on maps ( J z cos θ res , J z sin θ res ) $ (\sqrt{J_z}\,\cos{\theta_{\mathrm{res}}},\,\sqrt{J_z}\,\sin{\theta_{\mathrm{res}}}) $ we can see that flat BAN up and vCIR+ orbits have a common origin and circulate around the stable point that is offset from the origin. It should be added that in the absence of disturbances associated with the angle θres, these maps show only contours concentrated around the origin. When there are strong perturbations of ∼cos θres (t = 160 − 200), the maps show additional contours around a point offset from the origin.

thumbnail Fig. 3.

Evolution of instantaneous phase portraits on the (θz − θR, Jz) plane [ − π, π]×[0, 0.3] and on the J z ( cos ( θ z θ R ) , sin ( θ z θ R ) ) $ \sqrt{J_z}\left(\cos{(\theta_z-\theta_{\mathrm{R}}) },\sin{(\theta_z-\theta_{\mathrm{R}}})\right) $ plane [ − 0.5, 0.5]×[−0.5, 0.5] for three different HJ values. See Appendix B for a detailed description. The red box contains the phase portraits, which are discussed separately.

In a self-consistent way, the transformation of orbits in a distorted phase space transforms the phase space itself, which in turn changes the transformation paths of individual orbits. Figure 4 shows possible paths of orbit evolution under such a transformation. Orbits move from vCIR+ to BAN up or down or to vCIR− through vPAS up or down; from non-flat BAN up or down to vCIR−; and from flat BAN up to vCIR+, non-flat BAN up, or vCIR−, depending on a phase space structure. For example, at t = 160,  180 for HJ = −2.2, as seen in Figs. 3 and 4, there is no libration area near the stable point (θz − θR = 0, Jz ≈ 0.1), but there is a large area near (θz − θR = π, Jz ≈ 0.1). Thus, the orbits move from vCIR+ to BAN down or vCIR− through vPAS down, or to flat BAN up. Later, when a large number of vCIR+ orbits have transformed, conditions are created for a second stable point at θres = 0 (later t = 200). From this moment, the transformation of orbits begins to occur through both portals (θres = 0 and π), and the phase space gradually symmetrizes (Fig. 3, t = 300,  400).

thumbnail Fig. 4.

Phase portraits on the (θz − θR, Jz) plane for orbits along the major axis of the bar with fixed HJ = −2.2. These are the same portraits contained in the red box in Fig. 3. Portraits are numerically obtained at three points in time (t = 160,  180,  200). The bold arrows show the possible and forbidden (for a given HJ) paths of transformation of vCIR+ orbits (two left panels) and vCIR+ with BAN orbits (right panel) during the evolution of the phase space.

6. Discussion

The term buckling is employed to denote not only the rapid protrusion of the bar from the midplane, but also the asymmetric protrusion in the vertical direction (Raha et al. 1991; Martinez-Valpuesta et al. 2006; Sellwood & Gerhard 2020). This is believed to result in the formation of BP bulges. However, under certain conditions, BP bulges can also form completely symmetrically (Smirnov & Sotnikova 2018, 2019; Sellwood & Gerhard 2020; McClure et al. 2025). This difference in the morphology of the process forces researchers to look for different mechanisms that can facilitate the exit of the bar orbits from the midplane in both symmetric (resonant heating or resonant trapping) and asymmetric (instability) cases (Sellwood & Gerhard 2020).

In this Letter we want to emphasize that after the bar emergence and during its further evolution, the shape of the phase space changes due to perturbations of the Hamiltonian of the system. We leave aside the reason for these perturbations, but when the Hamiltonian of the system changes, both resonance trapping and heating can occur (Quillen 2006). This in turn leads to a redistribution of orbits in phase space. The system responds by further changing the Hamiltonian. Transformation of orbits and phase space occurs simultaneously and in a self-consistent manner.

Our simulations show that during buckling the asymmetric exit of orbits from the midplane is caused by perturbations of the Hamiltonian, both first- and second-order (Quillen et al. 2014), with both heating and trapping occurring simultaneously. If the amplitude of first-order Hamiltonian perturbations is large, the system admits only one fixed point in a wide range of values of the Jacobi integral (Fig. 3), and the exit of orbits from the midplane are asymmetric. Moreover, judging by the spike in the fraction of vPAS orbits in the first plot of Fig. 2, there is rather resonant heating during buckling. When first-order perturbations are absent or small (as in the case of a central concentration of mass, Smirnov & Sotnikova 2019; McClure et al. 2025), the transformation of the orbits occur through two fixed points – that is, almost symmetrically – and the resonant capture dominates. Zozulia et al. (2024b) show that at this stage the orbits can remain in vertical resonance for a long time.

Although Li et al. (2023), Łokas (2025) consider buckling to be a resonant phenomenon, we cannot agree with their interpretation of the orbital transformation. In their picture, the initial bisymmetric distortion of the Laplace plane at the ends of the bar (upward bending) leads to forced oscillations of the orbits with a vertical frequency of ν ≈ 2(Ω − Ωp). As a result, an excess of BAN up orbits is formed. Later, the distortion winds up and causes the transformation of banana-like orbits into pretzel-like one and gives rise to BAN down orbits.

In our picture, with the same initial distortion of the Laplace plane (due to the excess of vCIR+ up orbits), the role of a portal for transferring orbits to other regions of the phase space is played by a fixed point π, around which the BAN down orbits are concentrated from the very beginning. We consider that the transformation of orbits near this point occurs in the same way as in the second-order Hamiltonian model, with two fixed points (Quillen et al. 2014). Figure 3 shows that the height (Jz) of the fixed point π for the same value of the Jacobi integral decreases over time as the bar slows down, the Hamiltonian of the system changes, and the orbits are redistributed. The dropping down of fixed points reflects the phenomenon of separatrix shrinking (Quillen et al. 2014). In this case, the librating BAN orbits, which were initially captured from the vCIR+ region, can cross the separatrix, transforming into vCIR− orbits with high Jz values. We proceed from considerations of preserving the phase volume (Liouville’s theorem). As a result, if the orbits leave the circulation mode (vCIR+), the phase-space volume of the vCIR+ area decreases, and the volume of the areas where the orbits transfer increases. After a significant reduction in the area occupied by vCIR+ orbits, the conditions for the appearance of a second fixed point near zero arise.

We would also like to note that our classification of orbits does not contradict other classifications based, for example, on the ratio of vertical to radial oscillation frequencies, ωz/κ. The value of ωz/κ is often used to distinguish orbits from which different parts of the BP bulge are assembled (e.g., Portail et al. 2015; Parul et al. 2020). Our classification, based on the resonance angle behavior, paints a picture of the phase space in broad strokes. In it, for example, all orbits with a low frequency ratio fall into the region behind the separatrix, i.e., the region of vCIR− orbits, although among these there may be beautiful periodic orbits, such as “brezels” (Portail et al. 2015).

7. Conclusions

Zozulia et al. (2024b) studied the long-term and gradual evolution of the entire ensemble of orbits in the bar. Here, we focus on the early stages of a bar orbit evolution in the vertical direction, when all changes occur rapidly (buckling). Our analysis shows that differences in the behavior of the resonant angle correspond to different orbits in a phase space, which, in turn, create different structures in the xyz-space. During buckling, the perturbed phase space (the Hamiltonian of the system) evolves, and orbits change their types. The transformation of orbits leads to a transformation of phase space, which in turn leads to a redistribution of orbits with different behaviors in the resonant angle, θres. This happens in a self-consistent manner. As a result, the vertical shape of the bar changes.

Our main conclusions about what occurs during buckling are as follows:

1. The phase space of the bar is arranged as a layer cake. The lower layer, with small Jz (flat orbits), consists of circulating orbits with an increasing resonant angle (vCIR+). Next is the intermediate layer (the “filling”), containing banana-shaped orbits (BAN up and BAN down) that librate around 0 or π. Finally, the top layer contains orbits with large Jz that reduce the resonant angle (vCIR−). During the vertical growth of the bar, the shape of these layers changes, and the orbits are redistributed between them. An orbit can move from the lowest layer (vCIR+) to the middle layer (BAN) via resonant capturing, stay there for a while, and then move to the uppermost layer (vCIR−), or it can pass directly to vCIR− (vPAS) via resonant heating. In the case of a quiet evolution, the phase space (θres, Jz) of the bar is symmetric and changes slowly. Orbits move smoothly from one layer to another (Zozulia et al. 2024b) as the fixed points in phase space are lowered. The libration orbits either fall with them or gradually move to the upper layer behind the separatrix.

2. During buckling, the phase space (θres, Jz) is asymmetric, which means that the transfer of orbits to the higher layers is also asymmetric (Figs. 3 and 4).

3. At the beginning of buckling, there are practically no librating orbits near θres = 0 (BAN down), and the orbits move from vCIR+ to BAN up or vCIR− around θres = π (vPAS down). A new lowest layer of flat, librating orbits appears near θz − θR = 0. This layer is associated with the curvature of the minimum potential surface (Laplace plane) and exists as long as this curvature exists. Subsequently, a libration region forms at θres = 0 and high Jz. Orbits begin to be trapped in it, and a new transfer portal for vCIR+ orbits opens around θres = 0 (vPAS up). Gradually, during this transfer, the vCIR+ layer becomes thinner, the flat librating orbits disappear, and the phase space becomes symmetric. At this stage, the replenishment of the vCIR+ layer occurs at new orbits join the bar while its amplitude grows.


1

θres = π in our model.

2

Quillen et al. (2014), considering the general scheme of the phase space ( ( J z cos ϕ , J z sin ϕ ) $ (\sqrt{J_z}\,\cos{\phi},\,\sqrt{J_z}\,\sin{\phi}) $) transformation, employ another designation for the resonant angle ϕ = θz − 2(θφ − Ωpt). However, for the inner Lindblad resonance (ILR, bar orbits) 2θφ − Ωpt = 2θϕ ≈ θR, which means approximate equality θres = θz − θR ≈ ϕ.

3

When reducing a Hamiltonian system to a resonant one, averaging over fast angles occurs. Our procedure for calculating the action and angles already includes averaging over radial θR and vertical θz angles.

Acknowledgments

We acknowledge financial support from the Russian Science Foundation, grant no. 24-22-00376. We are also deeply grateful to an anonymous reviewer whose comments not only contributed to a much improved presentation of our results, but also generated a wide range of ideas for future research. We acknowledge the use of the AGAMA (Vasiliev 2019) and mpsplines (Ruiz-Arias 2022) python packages, without which the present work would not be possible.

References

  1. Binney, J. 2020, MNRAS, 495, 886 [NASA ADS] [CrossRef] [Google Scholar]
  2. Collier, A. 2020, MNRAS, 492, 2241 [Google Scholar]
  3. Combes, F., & Sanders, R. H. 1981, A&A, 96, 164 [NASA ADS] [Google Scholar]
  4. Combes, F., Debbasch, F., Friedli, D., & Pfenniger, D. 1990, A&A, 233, 82 [NASA ADS] [Google Scholar]
  5. Friedli, D., & Pfenniger, D. 1990, Eur. S. Obs. Conf. Workshop Proc., 35, 265 [Google Scholar]
  6. Kataria, S. K. 2024, MNRAS, 534, 3565 [Google Scholar]
  7. Li, X., Shlosman, I., Pfenniger, D., & Heller, C. 2023, MNRAS, 520, 1243 [NASA ADS] [CrossRef] [Google Scholar]
  8. Lichtenberg, A., & Lieberman, M. 1992, Regular and Chaotic Dynamics (New York: Springer-Verlag) [CrossRef] [Google Scholar]
  9. Łokas, E. L. 2019, A&A, 629, A52 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  10. Łokas, E. L. 2025, A&A, accepted [Google Scholar]
  11. Martinez-Valpuesta, I., Shlosman, I., & Heller, C. 2006, ApJ, 637, 214 [NASA ADS] [CrossRef] [Google Scholar]
  12. McClure, R. L., Beane, A., D’Onghia, E., Filion, C., & Daniel, K. J. 2025, MNRAS, 537, 1475 [Google Scholar]
  13. Navarro, J. F., Frenk, C. S., & White, S. D. M. 1996, ApJ, 462, 563 [Google Scholar]
  14. Parul, H. D., Smirnov, A. A., & Sotnikova, N. Y. 2020, ApJ, 895, 12 [Google Scholar]
  15. Pfenniger, D., & Friedli, D. 1991, A&A, 252, 75 [NASA ADS] [Google Scholar]
  16. Portail, M., Wegg, C., & Gerhard, O. 2015, MNRAS, 450, L66 [NASA ADS] [CrossRef] [Google Scholar]
  17. Quillen, A. C. 2002, AJ, 124, 722 [NASA ADS] [CrossRef] [Google Scholar]
  18. Quillen, A. C. 2006, MNRAS, 365, 1367 [Google Scholar]
  19. Quillen, A. C., Minchev, I., Sharma, S., Qin, Y.-J., & Di Matteo, P. 2014, MNRAS, 437, 1284 [Google Scholar]
  20. Raha, N., Sellwood, J. A., James, R. A., & Kahn, F. D. 1991, Nature, 352, 411 [Google Scholar]
  21. Ruiz-Arias, J. A. 2022, Solar Energy, 248, 121 [Google Scholar]
  22. Sellwood, J. A., & Gerhard, O. 2020, MNRAS, 495, 3175 [Google Scholar]
  23. Smirnov, A. A., & Sotnikova, N. Y. 2018, MNRAS, 481, 4058 [Google Scholar]
  24. Smirnov, A. A., & Sotnikova, N. Y. 2019, MNRAS, 485, 1900 [NASA ADS] [CrossRef] [Google Scholar]
  25. Sotnikova, N. Y., & Rodionov, S. A. 2003, Astron. Lett., 29, 321 [Google Scholar]
  26. Valencia-Enríquez, D., Puerari, I., & Chaves-Velasquez, L. 2023, MNRAS, 525, 3162 [Google Scholar]
  27. Vasiliev, E. 2019, MNRAS, 482, 1525 [Google Scholar]
  28. Zozulia, V. D., Smirnov, A. A., & Sotnikova, N. Y. 2024a, MNRAS, 529, 4405 [CrossRef] [Google Scholar]
  29. Zozulia, V. D., Smirnov, A. A., Sotnikova, N. Y., & Marchuk, A. A. 2024b, A&A, 692, A145 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

Appendix A: Orbital types

In a pendulum model, three modes of its behavior can be distinguished depending on the angle ϕ from the vertical, i.e., libration around a stable point and circulation with a decrease/increase in the angle (Lichtenberg & Lieberman 1992). These modes localize different regions of the phase space. The perturbed Hamiltonian dynamics near a resonance3 can be reduced to a first approximation to the pendulum model, where ϕ is equal to the resonant angle θres (Lichtenberg & Lieberman 1992; Binney 2020). By identifying the three modes of θres behavior, one can determine in which region of the phase space the orbit is located, at least for regular or sticky chaotic orbits that can change their behavior. It should be noted that, unlike classical mechanics systems, N-body models evolve by changing the potential (sometimes quite rapidly) and, therefore, the phase space. This means that orbits can change their behavior relative to the resonant angle. For this reason, we introduce an intermediate type of orbit behavior — “passage” — for orbits that pass from one circulation to another without getting stuck in resonance.

Obviously, the above reasoning applies to any resonant angle θres. The next step is to address the technical aspects of the matter. How can we identify the specific orbital types? To begin with, we must at least roughly understand the properties of the system: the position of the stable points and the boundaries within which libration regions may be located. In this Letter we consider two resonant angles 2θϕ − θR and θz − θR. The stable points of the first of them are 0 for x1 orbits (ILR in our notation) and π for x2 orbits, the libration region can be within ±π relative to these stable points. Similarly, the stable points of the second resonant angle are 0 or π. The exact resonances correspond to banana-shaped orbits, the tips of which are directed “up” and “down”, respectively.

It follows that the orbital type can be determined by the number of times the resonant angle takes a value equal to the angle of a fixed point and by whether it extends beyond the libration region. Zozulia et al. (2024b) we used a similar approach, but in the current study, we improved the methodology. If the system has stable points at θres = 0 ± π and π ± π, then the identification of types is as follows.

1, 2. Circulation of the resonant angle θres with positive θ ˙ res > 0 $ \dot{\theta}_{\mathrm{res}} > 0 $ or negative θ ˙ res < 0 $ \dot{\theta}_{\mathrm{res}} < 0 $ angular velocity. In this case, the resonant angle θres successively takes values that are multiples of π times in ascending or descending order, respectively.

3. Libration around a stable point. θres takes the value 0 or π at least three times without going beyond the boundaries of ±π. Thus, the orbit makes at least one revolution around a stable point in phase space. The orbit enters the resonance when it has Δθ12 left to pass to the stable point, where Δθ12 is the maximum deviation of the resonant angle relative to the stable point between the first and second passages through it. In fact, this represents the amplitude of the libration angle. The moment of exit from resonance is determined similarly.

4. Passage through the resonance is determined in the same way as libration, but θres takes the value 0 or π only twice. Formally, the orbit enters resonance, but does not make a full revolution around a stable point in phase space and changes the direction of its circulation. The initial and final times of passage are defined similarly to the libration.

Applying these criteria to the resonant angle 2θϕ − θR, we distinguish orbits in inner Lindblad resonance (ILR, x1), orbits with negative angular velocity Ω − κ/2 < Ωp, and all other orbits (x2, Residuals). We also divide the orbits in circulation into central and disk ones by the minimum distribution of the z-component of angular momentum Lz (Lz0).

In this Letter, ILR orbits are divided by the behavior of the resonant angle θz − θR and are briefly described in Sect. 3. The first and second types correspond to circulations with positive and negative angular velocity (vCIR+ and vCIR−), the third and fourth types correspond to librations (BAN) and passage (vPAS) near a stable point.

Appendix B: Sketches of phase portraits

The shape of phase space and the distribution of orbits within it completely determine the structure of a galaxy. However, the main challenge in studying galaxies in this way arises from the 6-dimensionality of this space. For a qualitative understanding of the features of the phase space, we need to reduce the dimensionality of the problem. Using the action-angle variables (JR, Jz, Lz, θR, θϕ, θz) helps us with this.

In Zozulia et al. (2024a,b) we describe the method for calculating the averaged action-angle variables that at least eliminates short-period variations of actions, frequencies, and angles. Thus, although not quite precisely (which, as follows from classical mechanics, requires a full Fourier transform), we eliminate the fast variables (angles θR and θz). This means that away from resonances we pass to three integrals of motion and a three-dimensional phase space. The problem becomes more complicated and interesting near the resonances.

We showed in Zozulia et al. (2024a) that using averaged action-angle variables (medium-term) allows us to reduce the number of dimensions to 4 near ILR and vILR. In this case, we assume that θ ˙ R , θ z ˙ θ z ˙ θ ˙ R , 2 θ ˙ ϕ θ ˙ R $ \dot{\theta}_{\mathrm{R}}, \dot{\theta_z} \gg \dot{\theta_z}-\dot{\theta}_{\mathrm{R}}, 2\dot{\theta}_\phi-\dot{\theta}_{\mathrm{R}} $, and then with a fixed value of the adiabatic invariant Jv = JR + Lz/2 + Jz or Jacoby integral HJ it is enough to consider the 4-dimensional space (Lz, θ1 = 2θϕ − θR, Jz, θ2 = θz − θR). In turn, if we consider orbits exactly along the bar, then we reduce the number of dimensions to two, implicitly reducing one of the equations of motion to the form θ 1 ˙ ( L z , J z , θ 1 = 0 , θ 2 ) = 0 $ \dot{\theta_1}(L_z, J_z, \theta_1 = 0, \theta_2) = 0 $. This means that, although qualitatively, it will allow us to explore the shape of phase space (θz − θR, Jz) and its evolution.

To investigate the phase space (θz − θR, Jz) we sketch the phase portraits. The algorithm for constructing them is divided into 2 steps. The first step is to search for orbits extended along the major axis of the bar. To do this, we solve a one-dimensional optimization problem. We launch orbits in a frozen rotating potential with coordinates (x, y = 0, z = z0, vx = 0, vy = HJ(x, z0),vz = 0), which correspond to (θR = 0, θz = 0, θϕ = 0). Here, z = z0 and the Jacobi integral HJ are fixed. We need to find the value of x such that the apocenters of the orbits lie on the major axis of the bar. To do this, we trace the evolution of the orbit for a given time (200 units, usually), and minimize the sum of the squares of the y-coordinate at the apocenter ∑y2(tapo) → min. The second step is to calculate the averaged medium-term θz, θR and Jz of the resulting, and plot the phase portrait (θz − θR, Jz) of this orbit with different initial values of z0.

It should be noted that only a part of the orbits in the bar is aligned exactly along its major axis, while most of them librate around it. The behavior of these orbits will be complicated as their phase space is four-dimensional. However, we believe that the general pattern of orbital evolution should be preserved. We expect only to see the splitting of the separatrix and the formation of a chaotic layer near it, which can accelerate the transition of orbits between different types.

All Figures

thumbnail Fig. 1.

Dynamic decomposition of the bar (ILR) at the beginning of buckling (t = 160, two top row) and at the peak of buckling (t = 180, two bottom row). For each time moment, the five leftmost snapshots on the xy-plane, covering [ − 2.5, 2.5]×[−2.5, 2.5], and on the xz-plane, covering [ − 2.5, 2.5]×[−1.5, 1.5], are shown. From left to right: All bar orbits (ILR); orbits with increasing resonant angle, θres = θz − θR (vCIR+); orbits in the vILR (BAN); orbits passing through it (vPAS); and orbits with decreasing θres (vCIR−). The next five snapshots show the same orbits on the xz-plane, but separated into “up” (near θres = 0) and “down” (near θres = π) types, and into secular Jz < 0.02 (very flat) and Jz > 0.02 (thick) for BAN orbits. For more information on the classification of these types, see Sect. 3 and Appendix A.

In the text
thumbnail Fig. 2.

The evolution of the fraction between different types of orbits from t = 100 to t = 300. In all plots, the buckling stage (t = 150 − 200) is highlighted by vertical lines. From left to right: among orbits in ILR (the bar), among ILR orbits captured in vILR (BAN), among ILR orbits passing through vILR (heating), among ILR orbits with increasing resonant angle θz − θR (vCIR+, flat orbits) and among ILR orbits with decreasing resonant angle θz − θR (vCIR−). The different types of orbits in plots, starting with the second, are shown in the same colors used in the first plot. All types correspond to those shown in Fig. 1.

In the text
thumbnail Fig. 3.

Evolution of instantaneous phase portraits on the (θz − θR, Jz) plane [ − π, π]×[0, 0.3] and on the J z ( cos ( θ z θ R ) , sin ( θ z θ R ) ) $ \sqrt{J_z}\left(\cos{(\theta_z-\theta_{\mathrm{R}}) },\sin{(\theta_z-\theta_{\mathrm{R}}})\right) $ plane [ − 0.5, 0.5]×[−0.5, 0.5] for three different HJ values. See Appendix B for a detailed description. The red box contains the phase portraits, which are discussed separately.

In the text
thumbnail Fig. 4.

Phase portraits on the (θz − θR, Jz) plane for orbits along the major axis of the bar with fixed HJ = −2.2. These are the same portraits contained in the red box in Fig. 3. Portraits are numerically obtained at three points in time (t = 160,  180,  200). The bold arrows show the possible and forbidden (for a given HJ) paths of transformation of vCIR+ orbits (two left panels) and vCIR+ with BAN orbits (right panel) during the evolution of the phase space.

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.