We numerically examine the corotation of two parallel quantized vortices in a selfgravitating BoseEinstein condensate BEC employing the GrossPitaevskiiPoisson equations. The longrange gravitationally attractive interaction allows the BEC to selfconfine without the need for external potentials while the densitydependence of the gravitational potential induces intriguing behaviors in the quantized vortices. The aim of this study is to provide a clue for understanding the corotation of two quantized vortices under the influence of gravitational interactions. The corotation of two quantized vortices is coupled with collective modes of the BEC which markedly differs from the behavior observed in typical BECs confined by an external potential. The rotational period increases linearly with the initial position from the center of the BEC. This deviation from the quadratic increase observed in a uniform BEC suggests that the gravitational interaction exerts a drag effect on the rotating quantized vortices. The two closely positioned quantized vortices rotate along elliptical orbits with radial fluctuations. However when the quantized vortices are initially positioned beyond a critical radius comparable to their core sizes their trajectory transitions into an outward spiral implying the onset of effective dissipation. Our findings demonstrate that the radial fluctuations of the quantized vortex resonate with the quadrupole mode of the BEC giving rise to a dissipation mechanism.
['Ultracold systems and matter waves', 'Condensed Matter - Other Condensed Matter', 'Astrophysics - Astrophysics of Galaxies']
Corotation of two quantized vortices coupled with collective modes in selfgravitating BoseEinstein condensates
{'Corotation of two quantized vortices coupled with collective modes in self-gravitating Bose-Einstein condensates': 'Kenta Asakawa 1, ∗ and Makoto Tsubota 1, 2, †  \n1 Department of Physics, Osaka Metropolitan University, 3-3-138 Sugimoto, Sumiyoshi-Ku, Osaka 558-8585, Japan 2 Nambu Yoichiro Institute of Theoretical and Experimental Physics (NITEP), Osaka Metropolitan University, 3-3-138 Sugimoto, Sumiyoshi-Ku, Osaka 558-8585, Japan (Dated: September 13, 2024)  \nWe numerically examine the corotation of two parallel quantized vortices in a self-gravitating Bose-Einstein condensate (BEC) employing the Gross-Pitaevskii-Poisson equations. The long-range gravitationally attractive interaction allows the BEC to self-confine without the need for external potentials, while the density-dependence of the gravitational potential induces intriguing behaviors in the quantized vortices. The aim of this study is to provide a clue for understanding the corotation of two quantized vortices under the influence of gravitational interactions. The corotation of two quantized vortices is coupled with collective modes of the BEC, which markedly differs from the behavior observed in typical BECs confined by an external potential. The rotational period increases linearly with the initial position from the center of the BEC. This deviation from the quadratic increase observed in a uniform BEC suggests that the gravitational interaction exerts a drag effect on the rotating quantized vortices. The two closely positioned quantized vortices rotate along elliptical orbits with radial fluctuations. However, when the quantized vortices are initially positioned beyond a critical radius comparable to their core sizes, their trajectory transitions into an outward spiral, implying the onset of effective dissipation. Our findings demonstrate that the radial fluctuations of the quantized vortex resonate with the quadrupole mode of the BEC, giving rise to a dissipation mechanism.', 'I. INTRODUCTION': "A quantized vortex, characterized by quantized circulation, plays a crucial role in the distinctive dynamics of Bose-Einstein condensates (BECs), superfluid 4 He, 3 He, and other quantum fluids [1]. It appears as a topological defect of the macroscopic wavefunction that describes the ordered state of the quantum fluid. Extensive theoretical and experimental research on the dynamics of the quantized vortex has been conducted in low-temperature physics, contributing significantly to the current understanding of some characteristic phenomena, such as the vortex lattice formation [2, 3] and quantum turbulence [4, 5]. Particularly, the behavior of the quantized vortices in the atomic BEC has been theoretically explored using the Gross-Pitaevskii (GP) model based on the mean-field approximation [6, 7].  \nInter-particle interactions are intricately connected to the behavior of the BEC and the motions of quantized vortices; the properties of the quantum fluid are significantly influenced by variations in the strength, range, and anisotropy of these interactions. In a typical atomic BEC, particles repel each other in short ranges due to contact interactions. This repulsive short-range interaction can be modulated through the Feshbach resonance to become attractive [9], potentially leading to instabilities in the BEC and the quantized vortices within it [10, 11]. Furthermore, in a dipolar BEC, dipole-dipole  \ninteractions, which vary in magnitude and direction depending on their orientation between particles, act over a wide range [12, 13]. The properties such as the long-range interactions and the switching between attraction and repulsion depending on the orientation between particles induce transverse instabilities in the BEC and straight quantized vortex [14, 15], as well as the shear motion in the Tkachenko mode exhibited by the striped vortex lattice [16]. The latest study investigates the impact of the dipole-dipole interactions on the generation of quantized vortices and the development of quantum turbulence by stirring the quasi-two-dimensional dipolar BEC with a rotating obstacle [17].  \nWhile low-temperature physics focuses on the atomic BEC where the contact and dipole-dipole interactions operate between atoms, astrophysics has paid attention to gravitationally interacting BECs composed of the ultralight scalar particle [18-25]. This is because such particles or BECs can be a plausible candidate for dark matter which constitutes the unidentified and primary component of the matter in the present universe. A BEC affected by the gravitational interaction is defined as a self-gravitating BEC. It is described by the non-local GP model coupled with the Newtonian gravitational field generated by the BEC itself. Since the gravitational interaction acts isotropically as an attractive long-range interaction, a self-gravitating BEC is confined by its intrinsic gravitational potential and can attain equilibrium independently, without the necessity of an external potential, distinguishing it from conventional atomic BECs [26, 27].  \nNotably, several optical methods have been proposed  \nto simulate the quantum fluid properties of such a selfgravitating BEC in laboratory systems. These include inducing the 'gravitation-like' interatomic interaction using off-resonant laser [27], and analogies with laser light propagating through thermo-optical media [28].  \nThe quantum-fluid dynamics of self-gravitating BECs with quantized vortices are intriguing from both the perspectives of low-temperature physics and astrophysics.  \nFrom the perspective of low-temperature physics, there is significant interest in the alterations in the behavior exhibited by BECs and their quantized vortices under the influence of isotropic long-range gravitational interactions. Their dynamics become nontrivial due to the density dependence of its own gravitational potential. For instance, the oscillation frequency of the collective modes in a self-gravitating BEC in response to small fluctuations depends on the total mass of the BEC [29, 30], reflecting the density dependence of the gravitational potential; this is different from an atomic BEC trapped by external potentials [8, 10]. When quantized vortices are present in a self-gravitating BEC, the change in the density profile caused by the vortex core deforms the gravitational potential [31], and both the BEC and quantized vortices move within this deformed potential. The intricate nonlinearity characteristic of self-gravitating BECs is likewise observed in dipolar BECs with long-range interactions, despite the differences in the anisotropy of these interactions.  \nOn the other hand, quantum fluid phenomena driven by quantized vortices in a self-gravitating BEC could serve as a touchstone for confirming the possibility that the BEC of ultralight scalar particles constitutes dark matter, as inferred from observations of galaxies and galaxy clusters. Astrophysical interest in this issue has been longstanding, encompassing concepts such as the possibility of the formation of vortex lattices in the ultralight dark matter BEC [32], comparisons between the critical angular velocity of a BEC with a central vortex and the rotation speed of galaxies [33], and the deformation of configurations due to quantized vortices [34]. Recently, three-dimensional numerical studies have explored the stability of the self-gravitating BEC with a single central vortex [35, 36], the gravitational wave emission during head-on collisions of two self-gravitating BECs with quantized vortices [31], and quantum turbulence resulting from the merger of multiple self-gravitating BECs [37-39].  \nAs an initial step toward understanding the dynamics of quantized vortices influenced by gravitational interactions in BECs, this study undertakes a numerical investigation of the orbital rotation of two quantized likecharged vortices in a three-dimensional self-gravitating BEC. There are numerous studies on the behavior of a few vortices in quantum fluids [1, 40-42], and it is wellknown that two parallel, straight vortices corotate within the quantum fluid, following the flow fields generated by each other. In a uniform quantum fluid at zero temperature, the two quantized vortices rotate while maintaining  \ntheir initial distance, and the rotational period increases proportionally to the square of the distance. On the other hand, when dissipation is present to reduce the rotational energy, the two parallel vortices spiral out as they rotate. Their corotational behavior is expected to be similarly exhibited even in a self-gravitating BEC, although the gravitational interaction, or the gravitational potential, is anticipated to influence their rotational speed and trajectory. Actually, recent studies on the rotation of two quantized vortices in a dipolar BEC, which also exhibits the long-range interaction, have revealed the non-circular rotational trajectories and the dependence of the rotational speed on the dipole-dipole interaction [43, 44].  \nWe compose this paper as follows. First, the formalism of the self-gravitating BEC and the numerical method based on it are shown in Sec.II. Our numerical results are shown in Sec.III. Finally, we conclude this paper in Sec.IV.", 'A. FORMULATION OF THE SELF-GRAVITATING BEC AND THE QUANTIZED VORTEX': 'When the self-gravitating BEC is composed of the Bose particles of the mass m and the s-wave scattering length a near zero temperature, its time evolution is governed by the Gross-Pitaevskii-Poisson (GPP) equations defined as  \ni ℏ ∂ψ ( r , t ) ∂t = -ℏ 2 2 m ∇ 2 ψ ( r , t ) + mV ( r , t ) ψ ( r , t ) + 4 π ℏ 2 a m | ψ ( r , t ) | 2 ψ ( r , t ) , (1)  \nand  \n∇ 2 V ( r , t ) = 4 πGm | ψ ( r , t ) | 2 , (2)  \nwhere, ψ ( r , t ) denotes the macroscopic wavefunction and V ( r , t ) represents the gravitational potential [2022]. The macroscopic wavefunction can be rewritten using the Madelung representation ψ ( r , t ) = √ ρ ( r , t ) /m exp[ iθ ( r , t )], where ρ ( r , t ) denotes the mass density and θ ( r , t ) is the phase. The velocity field is introduced as v ( r , t ) = ℏ ∇ θ ( r , t ) /m . Eq. (2) is the Poisson equation and determines the gravitational potential in Eq. (1) from the density profile of the BEC. Hence, the GPP equations exhibit that the self-gravitating BEC is constrained by and evolves according to its gravitational potential which depends on the density profile of the BEC itself.  \nAlthough the velocity field of the BEC is ordinarily irrotational, the quantized vortex defined as the topological deficit locally gives rise to the rotational flow in the BEC. For example, when the phase is given by θ ( r , t ) = tan -1 ( y/x ) using the Cartesian coordinates ( x, y, z ), the  \ncounter-clockwise velocity field is generated around the zaxis and its amplitude | v | ∝ 1 / √ x 2 + y 2 , which displays a straight vortex line along the z -axis. The quantized vortex is characterized by the discretized (quantized) circulation in units of κ = 2 π ℏ /m and the core size comparable to the coherence length ξ = √ m/ (8 πρ c a ), where ρ c denotes a characteristic value of the mass density. Note that the quantized vortex whose circulation exceeds the unit hardly appears because it is energetically unstable even in the self-gravitating BEC [35, 36].  \nThe total energy E ( t ) of the GPP equations (1) and (2) is given by the sum of the kinetic energy E k ( t ), gravitational energy E p ( t ), and contact interaction energy E i ( t ), namely E ( t ) = E k ( t ) + E p ( t ) + E i ( t ). These energy components are written by  \nE k ( t ) = ℏ 2 2 m ∫ d r |∇ ψ ( r , t ) | 2 , (3)  \nE p ( t ) = m 2 ∫ d r V ( r , t ) | ψ ( r , t ) | 2 , (4)  \nand  \nE i ( t ) = 2 π ℏ 2 a m ∫ d r | ψ ( r , t ) | 4 . (5)  \nUsing the Madelung representation, the kinetic energy is divided into two components: the quantum pressure component E q ( t ) and the velocity component E v ( t ) given by  \nE q ( t ) = ℏ 2 2 m 2 ∫ d r ∣ ∣ ∣ ∇ ( √ ρ ( r , t ) )∣ ∣ ∣ 2 , (6)  \nand  \nE v ( t ) = 1 2 ∫ d r ρ ( r , t ) v ( r , t ) 2 . (7)  \nE v is, furthermore, divided into the compressible part  \nE vc ( t ) = 1 2 ∫ d r ρ ( r , t ) v c ( r , t ) 2 , (8)  \nand incompressible part  \nE vi ( t ) = 1 2 ∫ d r ρ ( r , t ) v i ( r , t ) 2 , (9)  \nbased on the Helmholtz decomposition [45, 46]. Here, the compressible velocity field v c ( r , t ) and the incompressible one v i ( r , t ) are respetively satisfied that ∇× v c ( r , t ) = 0 and ∇· v i ( r , t ) = 0, and v c ( r , t ) + v i ( r , t ) = v ( r , t ). Although the deformation/oscillation of the BEC can contribute to both, the quantized vortices contribute only to the incompressible part.', 'B. NUMERICAL METHOD': "To carry out the numerical study of the self-gravitating BEC, we solve the dimensionless GPP equations,  \ni ∂ ˜ ψ (˜ r , ˜ t ) ∂ ˜ t = -1 2 ˜ ∇ 2 ˜ ψ (˜ r , ˜ t ) + ˜ V (˜ r , ˜ t ) ˜ ψ (˜ r , ˜ t ) +˜ a | ˜ ψ (˜ r , ˜ t ) | 2 ˜ ψ (˜ r , ˜ t ) , (10)  \nand  \n˜ ∇ 2 ˜ V (˜ r , ˜ t ) = | ˜ ψ (˜ r , ˜ t ) | 2 (11)  \nin a three-dimensional numerical box. Here, the dimensionless rules are respectively ˜ r = ˜ λ ( mc/ ℏ ) r , ˜ t = ˜ λ 2 ( mc 2 / ℏ ) t , ˜ ψ = { √ 4 πG ℏ / ( √ mc 2 ˜ λ 2 ) } ψ , ˜ V = V/ ( ˜ λ 2 c 2 ), and ˜ a = { ˜ λ 2 c 2 / ( mG ) } a , where c is the velocity of light in the vacuum [30, 47]. The criterion of these rules is the Compton wavelength, and the scaling parameter ˜ λ makes such dimensionless variables suitable without changing the formulations of Eqs. (10) and (11).  \nThe dimensionless GP equation of Eq. (10) is solved by the pseudo-spectrum method and the fourthorder Runge-Kutta method. Using the pseudo-spectrum method, the periodic boundary conditions must be applied to Eq. (10). On the other hand, the dimensionless Poisson equation of Eq. (11) can be solved by the Poisson solver in the Intel MKL. Considering the anisotropic deformation of the BEC near the center of the numerical box, we can employ the multipole expansion solution of the gravitational potential until the second-order terms, written by  \n˜ V (˜ r , ˜ t ) = -˜ M 4 π ˜ r -1 4 π ˜ r 3 ∑ i ∈{ x,y,z } ˜ P i ˜ r i -3 8 π ˜ r 5 ∑ i,j ∈{ x,y,z } ˜ I ij ˜ r i ˜ r j , (12)  \nas the boundary conditions for Eq. (11), where ˜ r i is i -component of ˜ r , e.g. ˜ r x = ˜ x , and ˜ r = | ˜ r | . ˜ M = ∫ d ˜ r | ˜ ψ (˜ r , ˜ t ) | 2 represents the dimensionless total mass, ˜ P i = ∫ d ˜ r ˜ r i | ˜ ψ (˜ r , ˜ t ) | 2 represents the i -component of the dimensionless dipole moment, and ˜ I i,j = ∫ d ˜ r (˜ r i ˜ r j -δ i,j ˜ r 2 / 3) | ˜ ψ (˜ r , ˜ t ) | 2 represents ( i, j )-component of the dimensionless quadrupole moment.  \nThe imaginary-time evolution of the GPP equations is applied to prepare the initial state. In this scheme, we initially set a suitable configuration and renew it by solving  \n-∂ ˜ ψ (˜ r , ˜ t ) ∂ ˜ t = -1 2 ˜ ∇ 2 ˜ ψ (˜ r , ˜ t ) + ˜ V (˜ r , ˜ t ) ˜ ψ (˜ r , ˜ t ) +˜ a | ˜ ψ (˜ r , ˜ t ) | 2 ˜ ψ (˜ r , ˜ t ) -˜ µ ( ˜ t ) ˜ ψ (˜ r , ˜ t ) , (13)  \nand the dimensionless Poisson equation of Eq. (11) until the system energically relaxes enough. The diffusion equation of Eq. (13) is obtained by a replacement  \n˜ t → -i ˜ t of Eq. (10) and the chemical potential ˜ µ ( ˜ t ) allows the system to relax while conserving the total mass. In the imaginary-time evolution, we employ the pseudo-spectrum method and Nesterov's accelerated gradient method [48-50] to implement that more rapidly.  \nWe prepare the initial state by positioning two parallel quantized vortices along the z-axis, aligned on the x-axis. The starting configuration of the imaginary-time evolution is written by ˜ ψ (˜ r , 0) = ˜ f (˜ r ) exp[ iθ (˜ r )]. The function ˜ f (˜ r ) is chosen such that it satisfies ˜ f (˜ x v, 1 , 0 , ˜ z ) = ˜ f (˜ x v, 2 , 0 , ˜ z ) = 0, where ˜ x v, 1 and ˜ x v, 2 are the x -coordinates of the quantized vortices, and the phase is set as  \nθ (˜ r ) = tan -1 ( ˜ y ˜ x -˜ x v, 1 ) +tan -1 ( ˜ y ˜ x -˜ x v, 2 ) . (14)  \nTwo conditions are imposed during the imaginary-time evolution of the GPP equations. First, the gravitational potential is symmetrized in xy , yz , and zx planes, namely ˜ V (˜ x, ˜ y, ˜ z, ˜ t ) = ˜ V (˜ x, ˜ y, -˜ z, ˜ t ), ˜ V (˜ x, ˜ y, ˜ z, ˜ t ) = ˜ V ( -˜ x, ˜ y, ˜ z, ˜ t ), and ˜ V (˜ x, ˜ y, ˜ z, ˜ t ) = ˜ V (˜ x, -˜ y, ˜ z, ˜ t ), in order to fix the center of mass of the BEC. Thus, we consider that the quantized vortices are symmetrically located in the yz plane, ˜ x v, 1 = -˜ x v, 2 = ˜ r o , where ˜ r o is the initial distance from the center of the BEC. Second, two initial identities,  \nRe [ ˜ ψ (˜ r , ˜ t ) ] = 0 when θ (˜ r ) = 0 , ± π, (15)  \nand  \nIm [ ˜ ψ (˜ r , ˜ t ) ] = 0 when θ (˜ r ) = ± π 2 (16)  \nare kept in order to fix the positions of the quantized vortices. These are the expansion of the conditions used in [51] which considers the case of one straight vortex. Adapting Eq. (14) and ˜ x v, 1 = -˜ x v, 2 = ˜ r o to Eqs. (15) and (16), these conditions in our situation are particularly rewritten as Re[ ˜ ψ (˜ r , ˜ t )] = 0 on a hyperbola ˜ x 2 -˜ y 2 = ˜ r 2 o and Im[ ˜ ψ (˜ x, 0 , ˜ z, ˜ t )] = 0.  \nFinally, we show the physical and numerical parameters for this work. We consider that the boson has the mass m = 3 × 10 -24 eV and the s-wave scattering length a ≈ 3 . 69 × 10 -81 m referred in [30, 35], and the total mass of the BEC is 3 × 10 12 M ⊙ which is comparable to the typical mass of the galaxy [52], where M ⊙ denotes the solar mass. The numerical box has the length 40 and the number of glids 256 3 , namely the spatial resolution is d ˜ x = 40 / 256 = 0 . 15625 or dx ≈ 0 . 10 kpc based on the dimensionless rule. The time resolution is 10 -3 . We determine the scaling parameter ˜ λ such that the dimensionless Thomas-Fermi radius ˜ R TF is 10. The Thomas-Fermi radius characterizes the size of the self-gravitating BEC and is written as R TF = π √ ℏ 2 a/ ( Gm 3 ) [20]. As a result, ˜ λ = 10 / ( πc ) √ Gm/a ≈ 3 . 30 × 10 -3 . The initial distance of the quantized vortex from the center of the BEC is set as r o /dx = 5 , 8 , 10 , 12 , 13 , 15 , 20 , 25 , 30 , 40, namely r o ≈ 0 . 5 , 0 . 8 , 1 . 0 , 1 . 2 , 1 . 3 , 1 . 5 , 2 . 0 , 2 . 5 , 3 . 0 , 4 . 0 in  \na unit of kilopersec. Note that all of r o is smaller than R TF ≈ 64 dx , which means that the vortices are put inner the BEC.", 'A. NON-CIRCULAR ORBITAL ROTATION OF TWO PARALLEL QUANTIZED VORTICES': 'Even in the self-gravitating BEC, two parallel quantized vortices co-rotate around the center of the BEC similarly to those in the conventional BEC as shown in Fig.1 and a corresponding video in the Supplemental Material [58]. It shows the mass density profiles in the plane z = 0 at various time when r o ≈ 2 . 0 kpc. Two dark-colored circles show the quantized vortices. Initially positioned on the x -axis as depicted in Fig. 1(a), they rotate at a constant angular velocity, completing a full revolution around 56 . 2 Myr, as illustrated in Fig. 1(e). Note that Fig. 1 also shows the deformation of the BEC, which is analyzed later.  \nThe rotation of two quantized vortices at a constant angular velocity is seen in the inset of Fig. 2. It is the time evolution of the rotational angle of two quantized vortices tan -1 ( y v, 1 /x v, 1 ) and tan -1 ( y v, 2 /x v, 2 ) when r o ≈ 1 . 0 kpc. These rotational angles linearly increase with the same tendency, ensuring both vortices rotate at the same constant angular velocity. Even though r o changes, the quantized vortices exhibit such a linear increase of the rotational angle.  \nAlthough the particular angular velocity or period characterizes the orbital rotation of two quantized vortices in the self-gravitating BEC similar to that in the usual BEC, its r o -dependence is different as shown in Fig. 2. In this figure, both red triangles and green crosses are not put on a gray dashed line, which shows the rotational period T rot does not quadratically increase with r o . Particularly, in r o ≲ 1 . 5 kpc they are put on a cyan straight line, and namely T rot ∝ r o . This differs from the rotational period of two vortices in a conventional BEC, namely T rot = (8 π 2 /κ ) r 2 o of Eq. (A6) which is derived in Appendix A, shown by the black dashed line. The extension of the period suggests the presence of an additional drag force acting on the vortices in the self-gravitating BEC. We infer that this drag force originates from gravitationally long-range interactions, analogous to the acceleration of the rotation of quantized vortices in a dipolar BEC [43, 44]. On the other hand, when r o ≳ 2 . 0 kpc, the rotational period does not increase significantly even if r o increases. We predict that it is caused by the finite size effect analogous to two vortices in a uniform BEC in a cylindrical container as shown in Appendix A.  \nIn the self-gravitating BEC, two quantized vortices display the bending of the vortex cores with corotating. Figure 3 and a corresponding video in the Supplemental Material [59] show the time evolution of two vortex cores when r o ≈ 1 . 0 kpc in the three-dimensional sys-  \n<latexit sha1\\_base64="Z/4dnOq/AkCehKTTpWGPvy5wLdY=">AAACZHichVHJSgNBEH0ZtxiXRIMgCBIMiqdQEc3iSfTi0USjgorMjJ04ZDZmJoEY/AG9Kh48KYiIn+HFH/CQHxDEo4IXD1YmERFRq+mu6tf1ql53K7auuR5RIyB1dHZ19wR7Q339A4PhyNDwumtVHFUUVEu3nE1FdoWumaLgaZ4uNm1HyIaiiw2lvNQ836gKx9Usc82r2WLHkEumVtRU2WMoV9uNxCmRzVImPRf7GSQT5FscbVuxItfYxh4sqKjAgIAJj2MdMlweW0iCYDO2gzpjDkeafy5wiBBzK5wlOENmtMxriXdbbdTkfbOm67NV7qLzdJgZwyQ90A290D3d0hO9/1qr7tdoaqmxV1pcYe+Gj0ZX3/5lGew97H+x/tTsoYiMr1Vj7baPNG+htvjVg7OX1fn8ZH2KLumZ9V9Qg+74Bmb1Vb3Kifw5QvwBn68c+z1Yn0kkU4lUbja+sNj+iiDGMIFpfu80FrCMFRS4r8AxTnAaeJT6pag00kqVAm1OFN9MGv8A2D2Kdg==</latexit>  \n<latexit sha1\\_base64="T74/OROw6Ku4+uEewTCQu8wCThE=">AAACdXichVHNLkNBFP56/ddPLzYSkTSKWNWpCCqRSGxsJC2KBGnuvaZ14/7l3mmTangAL2DBhgQRj2HjBSw8gliS2Fg4va2ICM5k5pz55nznfDOje5YZSKLHiNLU3NLa1t4R7ezq7ompvX3rgVvyDZEzXMv1N3UtEJbpiJw0pSU2PV9otm6JDX1/sXa+URZ+YLrOmqx4YsfWio5ZMA1NMpRXVTlPR9u2Jvd8u7pc8Q/zaoKS6fTMbGoy/jNIJSm0BBqWcdVrbGMXLgyUYEPAgeTYgoaAxxZSIHiM7aDKmM+RGZ4LHCLK3BJnCc7QGN3ntci7rQbq8L5WMwjZBnexePrMjGOUHuiGXuiebumJ3n+tVQ1r1LRU2Ot1rvDyseOB1bd/WTZ7ib0v1p+aJQqYDbWarN0LkdotjDq/fHDysjq3Mlodowt6Zv3n9Eh3fAOn/GpcZsXKKaL8AZ+vHP89WJ9MpqaT09mpxEK68RXtGMQwxvm9Z7CAJWSQ475lnOESV5E3ZUgZUcbqqUqkwenHN1MmPgCm7ZEp</latexit>  \n<latexit sha1\\_base64="ukh7u9QLbl17pgQPoA7D5OHVZfg=">AAACZHichVHLSsNAFD2Nr1ofrYogCFIsFVflVny0rkQ3LvuwWqhFkjhqaJqEJC3W4g/oVnHhSkFE/Aw3/oCL/oAgLhXcuPA2rYiIeoeZe+fMPfeemVEsXXNcooZP6ujs6u7x9wb6+gcGg6Gh4XXHrNiqyKmmbtp5RXaErhki52quLvKWLeSyoosNpbTSPN+oCtvRTGPNrVmiWJZ3DW1HU2WXofT+VihCsWSSEgtz4Z9BPEaeRdC2lBm6xia2YUJFBWUIGHA51iHD4VFAHASLsSLqjNkcad65wCECzK1wluAMmdESr7u8K7RRg/fNmo7HVrmLztNmZhhReqAbeqF7uqUnev+1Vt2r0dRSY6+0uMLaCh6NZd/+ZZXZu9j7Yv2p2cUOEp5WjbVbHtK8hdriVw/OXrKLmWh9ii7pmfVfUIPu+AZG9VW9SovMOQL8AZ+vHP49WJ+Jxedj8+nZyNJy+yv8GMckpvm9F7CEVaSQ474CxzjBqe9R6pdGpNFWquRrc0bwzaSJD9Y9inU=</latexit>  \n<!-- image -->  \n<!-- image -->  \n<!-- image -->  \n<!-- image -->  \nFIG. 1. The mass density profiles in the plane z = 0 when r o ≈ 2 . 0 kpc. (a) 0 Myr, (b) 14 . 0 Myr, (c) 28 . 1 Myr, (d) 42 . 1 Myr, and (e) 56 . 2 Myr. The box size is approximately 25 . 8 kpc.  \n<!-- image -->  \nFIG. 2. The relationship between the rotational period T rot and r o , compared to that of a conventional BEC given by Eq. (A6). The inset is the time evolution of the rotational angle of two quantized vortices in the plane z = 0 when r o ≈ 1 . 0 kpc.  \n<!-- image -->  \nFIG. 3. The vortex cores in the three-dimensional system at (a) 0 Myr, (b) 7 . 0 Myr, (c) 14 . 4 Myr, (d) 23 . 9 Myr, and (e) 30 . 6 Myr with their projections in xy , yz , and zx planes when r o ≈ 1 . 0 kpc. The red and green lines show the vortex cores.  \n<!-- image -->  \nFIG. 4. The trajectories of ( x v, 1 , y v, 1 ) (red line) and ( x v, 2 , y v, 2 ) (green line) shown by the orbital rotation of the quantized vortices. (a) r o ≈ 0 . 5 kpc, (b) r o ≈ 0 . 8 kpc, (c) r o ≈ 1 . 2 kpc, (d) r o ≈ 1 . 5 kpc, (e) r o ≈ 20 kpc, and (f) r o ≈ 2 . 5 kpc.  \n<!-- image -->  \ntem. Two vortex cores initially show straight lines along the z -axis in Fig. 3 (a). As the cores start to rotate, they undergo deformation, altering their bending direction as shown in Figs. 3 (b)-(e). Then, the vortex cores trace a non-circular trajectory as shown in Figs. 4 (a)(f). It shows the time evolution of the vortex positions in the z = 0 planes ( x v, 1 , y v, 1 ) and ( x v, 2 , y v, 2 ) as changing r o . The behavior can be divided into two types whether r o ≲ 1 . 0 kpc in Figs. 4 (a) and (b) or r o ≳ 1 . 2 kpc in Figs. 4 (c)-(f). In the former cases, the vortex cores follow elliptical orbits, whereas in the latter, they spiral outward with fluctuations.  \nFirst, let us examine the elliptical orbits in detail when r o is smaller than 1 . 0 kpc. We focus on two important terms; one is whether two orbits overlap or not, and the other is whether each orbit is closed or not. In Fig. 4 (a) corresponding to the case of r o ≈ 0 . 5 kpc, the elliptical  \n<latexit sha1\\_base64="hfhsQhk33vFWZkYgqVG1YglKUJg=">AAACeHichVHJSgNBEH0Z97hFvQhegsHtMlRicANB9OJFiEtUUAkzY6uDszHTCcSgH+APePAgCqLBz/DiD3jwE8SjgiAerEwiIqJW013Vr+tVve7WPcsMJNFDRKmrb2hsam6Jtra1d3TGurpXAzfvGyJruJbrr+taICzTEVlpSkuse77QbN0Sa/r+XOV8rSD8wHSdFVn0xJat7TrmjmlokqFcrEdOJ9MqHW3amtzz7dJC0T/MxRKkTk6mKD0a/xkkVQotgZpl3NgVNrENFwbysCHgQHJsQUPAYwNJEDzGtlBizOfIDM8FDhFlbp6zBGdojO7zusu7jRrq8L5SMwjZBnexePrMjGOA7qlMz3RHN/RI77/WKoU1KlqK7PUqV3i5zuPe5dd/WTZ7ib0v1p+aJXYwEWo1WbsXIpVbGFV+4eDkeXlqaaA0SBf0xPrP6YFu+QZO4cW4XBRLp4jyB3y+cvz3YDWlJsfUscV0Yma29hXN6EM/hvm9xzGDeWSQ5b5FnOEa5cibEleGlJFqqhKpcXrwzZTUBxcYkcg=</latexit>  \n<latexit sha1\\_base64="mr1TyqlhbxGvISrS9ui9bGZP55c=">AAACeHichVHJSgNBEH0Z97gk6kXwEgxul6ESQ1xAEL14EdyigkqYGVsdnI2ZTiAG/QB/wIMHURAVP8OLP+DBTxCPCoJ4sDKJiIhaTXdVv65X9bpb9ywzkEQPEaWuvqGxqbkl2trW3hGLd3atBG7BN0TOcC3XX9O1QFimI3LSlJZY83yh2bolVvW9mcr5alH4gek6y7LkiU1b23HMbdPQJEP5eLecTI+pqcMNW5O7vl2eK/kH+XiS1PHxNGVGEj+DlEqhJVGzeTd+iQ1swYWBAmwIOJAcW9AQ8FhHCgSPsU2UGfM5MsNzgQNEmVvgLMEZGqN7vO7wbr2GOryv1AxCtsFdLJ4+MxPop3u6pme6oxt6pPdfa5XDGhUtJfZ6lSu8fOyoZ+n1X5bNXmL3i/WnZoltjIVaTdbuhUjlFkaVX9w/fl6aWOwvD9A5PbH+M3qgW76BU3wxLhbE4gmi/AGfr5z4PVhJq6msml3IJKema1/RjF70YYjfexRTmMU8cty3hFNc4TrypiSUQWW4mqpEapxufDMl/QEjcZHO</latexit>  \n<latexit sha1\\_base64="Hcg1oAKlCIdy668YTx8c23wJz9Q=">AAACeHichVHJSgNBEH0Z97hFvQhegiEul1CJwQ0E0YsXwS0LqISZsdXB2ZjpBGLQD/AHPHgQBVHxM7z4Ax78BPGoIIgHK5OIiKjVdHf163rVr6s01zR8SfQQUhoam5pbWtvC7R2dXd2Rnt6s7xQ9XWR0x3S8vKb6wjRskZGGNEXe9YRqaabIaXvz1ftcSXi+4dhrsuyKTUvdsY1tQ1clQ4VIn5xJpxLJww1LlbueVVkseweFSIwSU1MpSo9FfzrJBAUWQ92WnMglNrAFBzqKsCBgQ7JvQoXPYx1JEFzGNlFhzGPPCO4FDhBmbpGjBEeojO7xusOn9Tpq87ma0w/YOr9i8vSYGUWc7umanumObuiR3n/NVQlyVLWUeddqXOEWuo/6V1//ZVm8S+x+sf7ULLGNyUCrwdrdAKn+Qq/xS/vHz6vTK/HKEJ3TE+s/owe65R/YpRf9YlmsnCDMDfiscvR3J8t9G0+ML6djs3P1VrRiAIMY4XpPYBYLWEKG3y3jFFe4Dr0pUWVYGa2FKqE6pw/fTEl9ABs3kco=</latexit>  \n<latexit sha1\\_base64="8nIEXY2NFQZ+nJjnYUSglljwERM=">AAACeHichVHJSgNBEH0Z97gk6kXwEgxul6ESNS4giF68CG5RQSXMjK0OzsZMJxCDfoA/4MGDKIiKn+HFH/DgJ4hHBUE8WJlEREStpruqX9eret2te5YZSKKHiFJTW1ff0NgUbW5pbYvF2ztWAjfvGyJruJbrr+laICzTEVlpSkuseb7QbN0Sq/reTPl8tSD8wHSdZVn0xKat7TjmtmlokqFcvFNOjmTU9OGGrcld3y7NFf2DXDxJ6vh4moaHEj+DlEqhJVG1eTd+iQ1swYWBPGwIOJAcW9AQ8FhHCgSPsU2UGPM5MsNzgQNEmZvnLMEZGqN7vO7wbr2KOrwv1wxCtsFdLJ4+MxPopXu6pme6oxt6pPdfa5XCGmUtRfZ6hSu8XOyoa+n1X5bNXmL3i/WnZoltjIVaTdbuhUj5FkaFX9g/fl6aWOwt9dE5PbH+M3qgW76BU3gxLhbE4gmi/AGfr5z4PVhJq6mMmlkYTk5NV7+iEd3owQC/9yimMIt5ZLlvEae4wnXkTUko/cpgJVWJVDmd+GZK+gMnkJHQ</latexit>  \nFIG. 5. (a) The time evolution of the distance of vortex 1 from the center of the BEC r v, 1 and the fitting results using the polar form of the ellipse equation. (b) The relationship between the initial distance of the vortex from the center of the BEC r o and the period of the oscillation estimated by fitting T rad , and the comparison to the rotational period of the vortex T rot and that of the quadrupole mode without the quantized vortex T q .  \n<!-- image -->  \ntrajectories just shift from each other. Conversely, as shown in Fig. 4 (b), for other cases such as r o ≳ 0 . 8 kpc, the two ellipses overlap. Additionally, as r o increases, the rotational orbits deviate from a perfectly elliptical shape, and as a result, the orbits are not strictly closed when r o ≳ 0 . 8 kpc. Such non-closed elliptical trajectories are also observed in a BEC with dipole-dipole interactions [43]. Hence, long-range interactions in a BEC may distort the trajectory of the corotating vortices from a closed circular form.  \nFigure 5 presents the quantitative analysis of the distance of each vortex from the center of the BEC r v,l = √ x 2 v,l + y 2 v,l , where l is 1 or 2, revealing the characteristics of elliptical orbits. As depicted in Fig. 5 (a), the quantized vortex exhibits radial oscillations. By fitting  \nthe data with the polar form of the ellipse equation, expressed as  \nr v,l = A (1 -ϵ 2 ) 1 -( -1) l ϵ cos ( 2 π T rad t ) , (17)  \nall results align with our numerical findings. Here, A represents the semi-major axis, ϵ denotes the eccentricity, and T rad is the period of the radial oscillation. This confirms quantitatively that the rotational orbits exhibit an elliptical form due to the radial oscillation of the quantized vortex. Additionally, the overlap of two ellipses and the non-closed trajectory can be explained by the eccentricity and the oscillatory period as follows. The eccentricity reaches a maximum of ϵ ≈ 0 . 09 when r o ≈ 0 . 5 kpc and decreases to ϵ ≈ 0 . 06 when r o ≳ 0 . 8 kpc, leading to the overlap of the orbits of co-rotating two vortices. According to Fig. 5 (b), which illustrates the r o -denpendence of T rad , it is comparable to the rotational period T rot when r o ≈ 0 . 5 kpc, but increases with r o , approaching a particular value. This increase in T rad relative to T rot as r o causes the rotational trajectory to become non-closed, due to phase discrepancies between the orbital rotation and the radial oscillation. Note that, as demonstrated later, the value toward which T rad converges with increasing r o is estimated to correspond to the period T q of the quadrupole mode in the self-gravitating BEC without vortices.  \nNext, we conduct a detailed analysis of the spiraling behavior of two quantized vortices when r o ≳ 1 . 2 kpc. Figure 4 (c) demonstrates that the rotational trajectories of the two quantized vortices exhibit a time-dependent transition from elliptical to spiral when r o ≈ 1 . 2 kpc; r v, 1 and r v, 2 initially oscillate until a certain point in time, after which they begin to increase while still oscillating. Such a time-dependent switching is absent when r o ≈ 1 . 5 kpc, as seen in Fig. 4 (d), where the two quantized vortices simply spiral outside while fluctuating. More intriguingly, for the cases where r o ≈ 2 . 0 kpc in Fig. 4 (e) and 2 . 5 kpc in Fig. 4 (f), two quantized vortices, after spiraling out sufficiently, then begin to spiral in.  \nThese spiral trajectories or increases of r v,l during rotation suggest the presence of dissipations acting on the quantized vortices, despite the fact that we consider the non-dissipative GPP model (see Eqs. (1) or (10)). We employ the vortex-point model with an effective dissipation analogous to mutual friction as a simplified approach to estimate the magnitude of these dissipative effects. For simplicity, we consider that two quantized vortices, with identical charges, are symmetrically positioned within a circular container of radius R . In this situation, the time evolution of the distance of each quantized vortex from the center of the BEC satisfies the following relation:  \nαt = Y [ r v,l ] (18)  \nas derived in Appendix A. Here, the parameter α quantitatively expresses the effective dissipation, and the func-  \nFIG. 6. The applications of the vortex-point model considering the two same-charged quantized vortices in a circular container to our numerical results. (a) The time evolution of a function Y [ r v,l ] which is explicitly written by Eq. (19). The solid lines show the numerical results obtained by substituting r v,l into Eq. (19), and the dotted lines demonstrate the linear fitting with those to estimate the mean value of α in Eq. (18). (b) The parameter α which shows the effective dissipation of the quantized vortex vs the initial distance of the quantized vortex from the center of the BEC, as indicated by black circles. The cyan line shows the coherence length ξ and the magenta line indicates twice that.  \n<!-- image -->  \nion Y [ r v,l ] is defined as  \nκ 2 π Y [ r v,l ] ≡ 4 √ 3 9 R 2 tan -1 ( √ 3 r 2 v,l R 2 ) -4 √ 3 9 R 2 tan -1 ( √ 3 r 2 o R 2 ) -r 2 v,l -r 2 o 3 . (19)  \nSubstituting the numerical results of r v,l into Eq. (19), the time evolution of Y [ r v,l ] can be obtained, as exhibited by solid curves in Fig. 6 (a). This figure allows us to categorize these solid curves into two groups; Y [ r v,l ] hardly increases when r o is small while it increases linearly, and dissipation works when r o is large. The mean increasing rate of Y [ r v,l ] corresponds to the magnitude of effective dissipation α , as indicated in Eq (18). By estimating the values of α through linear fitting as shown by the dotted lines in Fig. 6 (a), we observe the r o -dependence in Fig. 6 (b). The rapid increase between ξ and 2 ξ indicates that the critical distance at which the rotational trajectories transit from elliptical to spiral is comparable to the vortex core size. Therefore, we infer that inter-vortex interactions are closely related to the emergence of effective dissipation. It is noteworthy that effective dissipation also depends on time. For instance, when r o ≈ 2 . 0 kpc and 2 . 5 kpc, Y [ r v,l ] saturates around  \n<!-- image -->  \n/uni0028.9/uni0062.67 /uni0029.10  \nFIG. 7. The spectrum of r v, 1 using the FFT analyses when r o ≳ 1 . 2 kpc, namely the trajectory draws spiral. The solid lines respectively show the corresponding numerical results while the black dashed line is the frequency of the quadrupole mode in a self-gravitating BEC with no quantized vortex, represented as Eq. (20).  \n<!-- image -->  \nt = 40 Myr, and continues to fluctuate thereafter. Moreover, when r o ≈ 1 . 2 kpc and 1 . 3 kpc, the time-dependent transition of the trajectories from ellipses to spirals likely alters the magnitude of effective dissipation. However, in this study, we do not explore the time dependence of the effective dissipation in further detail.  \nWhile the spiral behavior of the two quantized vortices varies with r o , as shown in Figs. 4 (c)-(f), FFT analysis of the time evolution of r v,l reveals that their radial fluctuations share a common frequency. Figure 7 presents the spectra of r v, 1 when r o ≳ 1 . 2 kpc, obtained from the FFT analysis. These spectra consistently exhibit a peak corresponding to the quadrupole mode in a self-gravitating BEC without vortices, whose period is  \ngiven by  \nT q = 2 π √ 5 π ( π 2 -6) 3 √ ℏ 6 a 3 G 5 m 9 1 √ M (20)  \nassuming the total mass M is sufficiently large [29, 30]. Hence, r v,l oscillates with a frequency comparable to that of the quadrupole mode even though r o just changes in the situations where the quantized vortices spiral outward. We can conclude that, as r o increases, the period of the radial oscillation T rad becomes smaller than the rotational one T rot and approaches T q .  \nIn all the figures presented above, the numerical results for cases where r o ≈ 3 . 0 kpc and 4 . 0 kpc are not contained. This is due to the unexpected emergence of the ghost vortices near the edge of the BEC, which hinders the accurate determination of the coordinates of the two quantized vortices. Figures 8 (a)-(e) and a corresponding video in the Supplemental Material [60] depict the time evolution of the mass density profiles, sliced at the z = 0 plane, where r o ≈ 3 . 0 kpc. In the early stage demonstrated in Fig. 8 (b), surface fluctuations induced by the rotation of the quantized vortices cause partial concavities at the edge of the BEC, compared to the initial profile in Fig. 8 (a). These concavities evolve into two ghost vortices in Fig. 8 (c), and their number gradually increases to 4-6, as illustrated in Figs. 8 (d) and (e). Furthermore, according to Fig. 8 (e), two quantized vortices complete approximately one rotation at t ≈ 58 Myr when r o ≈ 3 . 0 kpc. This rotational period is shorter than that for r o ≈ 2 . 5 kpc, as shown in Fig. 2, supporting the probability that finite-size effects accelerate the rotational speed of the quantized vortices, even in a self-gravitating BEC.', 'B. DEFORMATION, ROTATION, AND OSCILLATION OF THE BEC': 'Figure 1 and Figure 8 also demonstrate that as the two quantized vortices corotate, the BEC deforms into an ellipsoidal shape and rotates around the z -axis. This ellipsoidal deformation and rotation are quantitatively analyzed by tracking the effective sizes in each direction X ( t ), Y ( t ), and Z ( t ) defined as  \nX ( t ) = √ m M ∫ d r x 2 | ψ ( r , t ) | 2 , (21)  \nY ( t ) = √ m M ∫ d r y 2 | ψ ( r , t ) | 2 , (22)  \nZ ( t ) = √ m M ∫ d r z 2 | ψ ( r , t ) | 2 . (23)  \nand  \nIn our numerical results, the time evolutions of these semi-axes exhibit complex oscillations, indicating that the superposition of collective modes drives the deformation and rotation within the BEC.  \nWe focus on three collective modes: the spherical breathing (monopole) mode, the quadrupole mode along the z -axis, and the quadrupole mode in the xy plane. These collective modes are characterized by specific quantities of the ellipsoidal shape: the mean radius ( XYZ ) 1 / 3 and two aspect ratios of the BEC Z/ √ X 2 + Y 2 and ( X 2 -Y 2 ) / ( X 2 + Y 2 ). As shown in Fig. 9, the time evolution of these quantities exhibits different oscillations. These figures reveal that each oscillation has a different dependence on r o . As r o increases, the oscillation of the mean radius results from the superposition of several modes, the oscillation of Z/ √ X 2 + Y 2 remains monotonic and independent of r o , and the oscillation of ( X 2 -Y 2 ) / ( X 2 + Y 2 ) converges to a specific mode. According to Fig. 9 (a), the size of the BEC decreases with increasing r o [61]. On the other hand, Figs. 9 (b), (c1), and (c2) demonstrate that as r o increases, the initial configuration becomes distorted from the spherical shape, leading to an amplification of the oscillations in both aspect ratios. The most intriguing aspect is that, as shown in Fig. 9 (c2), is the gradual increase in the amplitude of 1 . 2 kpc ≲ r o ≲ 2 . 0 kpc, despite the system being non-dissipative. This suggests that the quadrupole mode in the xy plane contributes to the effective dissipation of the corotating quantized vortex.  \nAs observed in Fig. 9, the presence of quantized vortices also affects the frequency of these collective modes, particularly the quadrupole mode in the xy plane. Figure 10 shows the maximum spectra peaks obtained by the FFT analysis of the oscillations in Fig. 9, illustrating the relationship between r o and the periods of each collective mode. The breathing mode shows that, while its period remains largely unaffected by r o , the presence of quantized vortices causes a shift in its frequency from that of the no-vortex case, which is given by  \nT b = 2 π √ 2 π ( π 2 -6) 3 √ ℏ 6 a 3 G 5 m 9 1 √ M (24)  \n[29, 30]. This increase in the period due to the presence of quantized vortices qualitatively aligns with the analytical results concerning a quantized vortex at the center of a self-gravitating BEC [53]. While the quadrupole mode along the z -axis remains unaffected by the quantized vortices, the quadrupole mode in the xy -plane exhibits an increase in its period with r o , approaching T q expressed in Eq. (20). We infer that the r o -dependence of the quadrupole mode in the xy plane is linked to the rotation of the BEC, caused by the splitting of two degenerate frequencies due to the presence of the quantized vortices. It is well-known that such rotation occurs when a quantized vortex is situated at the center of a trapped BEC, and its angular velocity is determined by the difference between the two shifted frequencies [8].  \n<latexit sha1\\_base64="Z/4dnOq/AkCehKTTpWGPvy5wLdY=">AAACZHichVHJSgNBEH0ZtxiXRIMgCBIMiqdQEc3iSfTi0USjgorMjJ04ZDZmJoEY/AG9Kh48KYiIn+HFH/CQHxDEo4IXD1YmERFRq+mu6tf1ql53K7auuR5RIyB1dHZ19wR7Q339A4PhyNDwumtVHFUUVEu3nE1FdoWumaLgaZ4uNm1HyIaiiw2lvNQ836gKx9Usc82r2WLHkEumVtRU2WMoV9uNxCmRzVImPRf7GSQT5FscbVuxItfYxh4sqKjAgIAJj2MdMlweW0iCYDO2gzpjDkeafy5wiBBzK5wlOENmtMxriXdbbdTkfbOm67NV7qLzdJgZwyQ90A290D3d0hO9/1qr7tdoaqmxV1pcYe+Gj0ZX3/5lGew97H+x/tTsoYiMr1Vj7baPNG+htvjVg7OX1fn8ZH2KLumZ9V9Qg+74Bmb1Vb3Kifw5QvwBn68c+z1Yn0kkU4lUbja+sNj+iiDGMIFpfu80FrCMFRS4r8AxTnAaeJT6pag00kqVAm1OFN9MGv8A2D2Kdg==</latexit>  \n<latexit sha1\\_base64="T74/OROw6Ku4+uEewTCQu8wCThE=">AAACdXichVHNLkNBFP56/ddPLzYSkTSKWNWpCCqRSGxsJC2KBGnuvaZ14/7l3mmTangAL2DBhgQRj2HjBSw8gliS2Fg4va2ICM5k5pz55nznfDOje5YZSKLHiNLU3NLa1t4R7ezq7ompvX3rgVvyDZEzXMv1N3UtEJbpiJw0pSU2PV9otm6JDX1/sXa+URZ+YLrOmqx4YsfWio5ZMA1NMpRXVTlPR9u2Jvd8u7pc8Q/zaoKS6fTMbGoy/jNIJSm0BBqWcdVrbGMXLgyUYEPAgeTYgoaAxxZSIHiM7aDKmM+RGZ4LHCLK3BJnCc7QGN3ntci7rQbq8L5WMwjZBnexePrMjGOUHuiGXuiebumJ3n+tVQ1r1LRU2Ot1rvDyseOB1bd/WTZ7ib0v1p+aJQqYDbWarN0LkdotjDq/fHDysjq3Mlodowt6Zv3n9Eh3fAOn/GpcZsXKKaL8AZ+vHP89WJ9MpqaT09mpxEK68RXtGMQwxvm9Z7CAJWSQ475lnOESV5E3ZUgZUcbqqUqkwenHN1MmPgCm7ZEp</latexit>  \n<latexit sha1\\_base64="ukh7u9QLbl17pgQPoA7D5OHVZfg=">AAACZHichVHLSsNAFD2Nr1ofrYogCFIsFVflVny0rkQ3LvuwWqhFkjhqaJqEJC3W4g/oVnHhSkFE/Aw3/oCL/oAgLhXcuPA2rYiIeoeZe+fMPfeemVEsXXNcooZP6ujs6u7x9wb6+gcGg6Gh4XXHrNiqyKmmbtp5RXaErhki52quLvKWLeSyoosNpbTSPN+oCtvRTGPNrVmiWJZ3DW1HU2WXofT+VihCsWSSEgtz4Z9BPEaeRdC2lBm6xia2YUJFBWUIGHA51iHD4VFAHASLsSLqjNkcad65wCECzK1wluAMmdESr7u8K7RRg/fNmo7HVrmLztNmZhhReqAbeqF7uqUnev+1Vt2r0dRSY6+0uMLaCh6NZd/+ZZXZu9j7Yv2p2cUOEp5WjbVbHtK8hdriVw/OXrKLmWh9ii7pmfVfUIPu+AZG9VW9SovMOQL8AZ+vHP49WJ+Jxedj8+nZyNJy+yv8GMckpvm9F7CEVaSQ474CxzjBqe9R6pdGpNFWquRrc0bwzaSJD9Y9inU=</latexit>  \nFigure 12 suggests that the effective dissipation of the quantized vortices is associated with its resonance with the quadrupole mode. As shown in Figs. 12 (a) and (b), the incompressible part E vi ( t ) defined as Eq. (9) predominantly governs the velocity-field contribution of the kinetic energy E v ( t ). However, in Figs. 12 (b), E vi ( t ) exhibits a slight decrease when 1 . 2 kpc ≲ r o ≲ 2 . 5 kpc where resonance between the quantized vortices and the quadrupole mode in the xy plane occurs. Since the system conserves the total energy, this reduction suggests that the energy of the quantized vortices is gradually redistributed into other energy components during their corotation. In contrast, the time evolution of  \n<!-- image -->  \n<!-- image -->  \n<!-- image -->  \nFIG. 8. The mass density profiles in the plane z = 0 when r o ≈ 3 . 0 kpc. (a) 0 Myr, (b) 19 . 2 Myr, (c) 38 . 3 Myr, (d) 44 . 7 Myr, and (e) 57 . 5 Myr. The box size is approximately 25 . 8 kpc.  \n<!-- image -->  \nFIG. 9. Time evolutions of quantities that characterize an ellipsoidal shape of the BEC when r o changes. (a) the average radius ( XYZ ) 1 / 3 , (b) the aspect of the BEC along the z -axis Z/ √ X 2 + Y 2 , (c1) and (c2) the aspect of the BEC in the xy planes ( X 2 -Y 2 ) / ( X 2 + Y 2 ).  \n<!-- image -->  \nAs mentioned above, the oscillations of the mean radius and the aspect ratio in the xy plane are superimposed across several modes depending on r o . In particular, the latter is likely attributable to the lifting of the degeneracy in the quadrupole mode. However, we were unable to detect these frequencies in the present study, probably because the lower split frequency does not lie within the time window covered by our numerical calculations. Further investigations will be pursued as part of our future work.', 'C. ENERGY TRANSFER TRIGGERED BY RESONANCE BETWEEN QUANTIZED VORTICES AND QUADRUPOLE MODE': 'When the two corotating quantized vortices trace spiral trajectories, their radial oscillations resonate with the quadrupole mode of the BEC. As shown in Fig. 11, T rad increases with r o and gradually approaches the period of ( X 2 -Y 2 ) / ( X 2 + Y 2 ). As a result, when r o ≳ 1 . 2 kpc, or two quantized vortices spiral outward, the radial oscillation of the quantized vortices has a similar period to the  \nFIG. 10. r o -dependence of the periods corresponding to a maximum peak of spectra of each oscillation (black circles: the mean radius, red triangles: the aspect ratio of the BEC along the z -axis, green squares: the aspect ratio of the BEC in the xy plane). The solid lines show the period of collective modes in a self-gravitating BEC without quantized vortices: magenta represents the breathing mode and cyan represents the quadrupole mode.  \n<!-- image -->  \noscillation of the aspect ratio ( X 2 -Y 2 ) / ( X 2 + Y 2 ), corresponding to the quadrupole mode in the xy plane. Thus, when two quantized vortices rotate in this system, their radial fluctuations couple with the quadrupole mode.  \n<latexit sha1\\_base64="TXzYdofTmePplHtx+6l+SYd+M1I=">AAACeHichVHJSgNBEH0Ztxi3qBfBSzC4XYaeINEIguDFi+AWDSQhzIydOGQ2ZjqBGOIH+AMePIiCqPgZXvwBD36CeFQQxIOVSURE1Gq6q/p1varX3ZprGr5g7CEkdXR2dfeEeyN9/QODQ9HhkR3fqXg6T+uO6XgZTfW5adg8LQxh8ozrcdXSTL6rlVea57tV7vmGY2+LmsvzllqyjaKhq4KgQnRULCkpOXGYs1Sx71n1tZrXKETjTE6l5heUROxnoMgssDjatu5EL5HDHhzoqMAChw1BsQkVPo0sFDC4hOVRJ8yjyAjOORqIELdCWZwyVELLtJZol22jNu2bNf2ArVMXk6ZHzBgm2T27Zs/sjt2wR/b+a616UKOppUZea3G5Wxg6Gtt6/ZdlkRfY/2L9qVmgiIVAq0Ha3QBp3kJv8asHx89bi5uT9Sl2zp5I/xl7YLd0A7v6ol9s8M0TROgDPl859nuwk5CVpJzcmIsvp9pfEcY4JjBD7z2PZaxiHWnqW8MprnAdepNi0rQ020qVQm3OKL6ZlPgAVwCR4Q==</latexit>  \n(c)  \n<latexit sha1\\_base64="1vTHlch8XNKmSKy3CaFEf3GhGaY=">AAACeHichVHJSgNBEH0Z97hFvQhegsHtMnQS0UQQBC9ehCxGBRPCzNjqkNmY6QRi0A/wBzx4EAXR4Gd48Qc85BPEo4IgHqxMIiKiVtNd1a/rVb3uVh1D9wRjjYDU0dnV3dPbF+wfGBwaDo2Mbnp22dV4TrMN291WFY8busVzQhcG33ZcrpiqwbfU0mrzfKvCXU+3rQ1RdXjBVPYtfU/XFEFQMTQmluMJOX6cNxVx4Jq19ap7VAxFmJxMLiaisfDPICoz3yJoW8oOXSOPXdjQUIYJDguCYgMKPBo7iILBIayAGmEuRbp/znGEIHHLlMUpQyG0ROs+7XbaqEX7Zk3PZ2vUxaDpEjOMKfbA6uyZ3bNb9sjef61V82s0tVTJqy0ud4rDJ+PZ139ZJnmBgy/Wn5oF9pDwteqk3fGR5i20Fr9yePqcXcpM1abZJXsi/Reswe7oBlblRbtK88wZgvQBn68c/j3YjMnRBXkhPR9ZSba/ohcTmMQsvfciVrCGFHLUt4pz3KAeeJPC0ow010qVAm3OGL6ZFPsAWx6R4w==</latexit>  \nt  \n= 38  \n.  \n3 Myr  \n<latexit sha1\\_base64="M927hv/DSRwkOey+pxQsIZwVdIc=">AAACeHichVFLSgNBEH0Z//GTaDaCm2Dwtxk6ISRGEAQ3bgQ1RgWVMDN2kiHzY6YTiCEewAu4cCEKouIx3HgBFx5BXEYQxIWVSURE1Gq6q/p1varX3apj6J5g7DEgdXX39Pb1DwQHh4ZHQuHRsS3Prrgaz2m2Ybs7quJxQ7d4TujC4DuOyxVTNfi2Wl5unW9XuevptrUpag7fN5WipRd0TREE5cMRsZhMyumjPVMRJdesr9bcRj4cY3Imk56PJ6I/g7jMfIuhY2t2+Ap7OIANDRWY4LAgKDagwKOxizgYHML2USfMpUj3zzkaCBK3QlmcMhRCy7QWabfbQS3at2p6PlujLgZNl5hRTLEHdsOa7J7dsif2/mutul+jpaVGXm1zuZMPHY9nX/9lmeQFSl+sPzULFDDva9VJu+MjrVtobX718KSZXdiYqk+zC/ZM+s/ZI7ujG1jVF+1ynW+cIkgf8PnK0d+DrYQcT8mp9WRsKdP5in5MYBKz9N5pLGEFa8hR3xrOcI2bwJsUlWakuXaqFOhwIvhmUuIDXSaR5A==</latexit>  \n<latexit sha1\\_base64="hlFcoFI0xaDIB4pvZTvEz6RlnAc=">AAACeHichVHJSgNBEH0Z97gk6kXwEgxul6ETjDGCIHjxIrglEUwIM2MnGTIbM51ADPoB/oAHD6IgKn6GF3/Ag58gHhUE8WBlEhERtZruqn5dr+p1t+oYuicYewhIHZ1d3T29fcH+gcGhUHh4JOPZVVfjac02bHdHVTxu6BZPC10YfMdxuWKqBs+qlZXmebbGXU+3rW1Rd3jeVEqWXtQ1RRBUCI+KpURSThzmTEWUXbOxVncPCuEok1Op5EIsHvkZxGTmWxRtW7fDl8hhDzY0VGGCw4Kg2IACj8YuYmBwCMujQZhLke6fcxwgSNwqZXHKUAit0Fqi3W4btWjfrOn5bI26GDRdYkYwye7ZNXtmd+yGPbL3X2s1/BpNLXXyaovLnULoaGzr9V+WSV6g/MX6U7NAEQu+Vp20Oz7SvIXW4tf2j5+3FjcnG1PsnD2R/jP2wG7pBlbtRbvY4JsnCNIHfL5y5PcgE5dj8/L8xlx0OdX+il6MYwIz9N5JLGMV60hT3zpOcYXrwJsUkaal2VaqFGhzRvHNpPgHYUmR5g==</latexit>  \nFIG. 11. The relationship between r o and three periods. Black circles show the period rotating two quantized vortices T rot , green squares show that oscillating in the radial direction T rad , and red triangles show that of ( X 2 -Y 2 ) / ( X 2 + Y 2 ), namely the period of the quadrupole mode in the xy plane.  \n<!-- image -->  \nFIG. 12. Time evolutions of kinetic energy components when r o changes. (a) the velocity-field component E v ( t ), (b) incompressible part E vi ( t ), (c1) and (c2) compressible part E vc ( t ).  \n<!-- image -->  \nthe compressible part E vc ( t ) in Figs. 12 (c1) and (c2) shows distinct qualitative changes depending on whether 1 . 2 kpc ≲ r o ≲ 2 . 5 kpc or not. The former case in Fig. 12 (c2) shows E vc ( t ) increases with fluctuating and saturates in late time, which implies that phonons are excited in the BEC until saturation. On the other hand, the latter one in Fig. 12 (c1) does not cause the rise of E vc ( t ). Hence, oscillations or deformations of the BEC are time-dependently amplified when the effective dissipation works in the quantized vortices, likely caused by the energy transfer by phonons from the incompressible part. Such time-dependent amplification of oscillation is also shown by ( X 2 -Y 2 ) / ( X 2 + Y 2 ), and therefore we  \nconclude that the quadrupole mode in the xy plane receives the energy from vortices through emitted phonons. Notably, the time evolution of the gravitational component of Eq. (4) and the contact-interaction component of Eq. (5) exhibit the complicated oscillation similar to the mean radius of the BEC as seen in Fig. 9 (a), which reflects that the configuration is composed of the competition between two self-interactions.  \nFigure 12 (a) is likely to indicate that the corotation of the two quantized vortices becomes coupled with the quadrupole mode in the xy plane when r o becomes sufficiently large. For cases where r o ≳ 3 . 0 kpc, the time evolution of E vc ( t ) exhibits more pronounced oscillations compared to when r o ≲ 2 . 5 kpc, yet the increase observed in Fig. 12 (b) does not manifest. We infer that the rotating quantized vortices near the BEC surface interfere with the quadrupole mode in the xy plane, making it challenging to discern whether the deformation and rotation are primarily driven by the quadrupole mode or the quantized vortices.', 'IV. CONCLUSIONS AND DISCUSSIONS': 'We numerically investigate the corotation of two parallel quantized vortices in a self-gravitating BEC using the three-dimensional GPP equation. The gravitational potential that confines the BEC is influenced by the motion of the BEC and the presence of quantized vortices through its density dependence. As a result, the dynamics of quantized vortices become complex and highly intriguing. This work seeks to elucidate the effects of gravitational interactions on the motion of quantized vortices within the BEC and to uncover non-trivial quantum hydrodynamical phenomena driven by self-gravity.  \nThe gravitational interaction exerts a drag force on the rotating quantized vortices, thereby prolonging their rotational period. When the two quantized vortices are initially positioned near the center of the BEC, their rotational period increases linearly with their initial distance from the center, differing from the quadratic increase observed in a uniform BEC. This difference suggests that the gravitational interaction acts as a drag force on the co-rotating vortices and reduces their rotational velocity. Conversely, when the vortices are initially positioned near the surface of the BEC, the rotational period does not increase significantly due to finite-size effects.  \nThe corotating quantized vortices in the selfgravitating BEC exhibit radial oscillations. As the initial distance of the vortices from the center increases, this oscillation period extends and converges towards that of the quadrupole mode in a self-gravitating BEC without quantized vortices. Compared to the rotational period, the period of radial oscillations aligns when the two vortices are sufficiently close but becomes relatively shorter as they are further separated. The resulting discrepancy between the two periods leads to a non-closed rotational trajectory of two vortices in the xy plane. The most in-  \ntriguing aspect is that the shape of the trajectory qualitatively changes depending on the initial position of the quantized vortices. When the vortices are initially positioned near the center, they rotate along an elliptical orbit; otherwise, they spiral outward with fluctuations. The critical distance at which the trajectory transitions from elliptical to spiral is comparable to the core size of the quantized vortices, suggesting that inter-vortex interactions play a role in this transition. Notably, the vortex cores bend from a straight line along the rotational axis, but we infer that this deformation does not significantly contribute to the trajectory transition based on an analogy with the corotation in a quasi-two-dimensional dipolar BEC [43, 44].  \nAs the two quantized vortices corotate, the selfgravitating BEC deforms and rotates due to the superposition of multiple collective modes. The effective sizes of the BEC in each direction allow us to identify three distinct collective modes: the breathing mode, the quadrupole mode along the z -axis, and the quadrupole mode in the xy -plane. Focusing on the periods of these collective modes, their dependencies on the initial position of the quantized vortices are different. The periods of the breathing mode and the quadrupole mode along the z -axis remain nearly constant, regardless of the initial vortex positions. In contrast, the farther the quantized vortices are initially located from the center, the longer the period of the quadrupole mode becomes, gradually approaching the period of the quadrupole mode observed in the absence of quantized vortices. The amplitudes of these collective modes are also influenced by the initial positions of the quantized vortices. Particularly, the quadrupole mode in the xy plane becomes progressively amplified over time when the corotation of the two quantized vortices results in outward spiral trajectories.  \nWe conclude that the effective dissipation of the quantized vortices is triggered by the coupling between the radial fluctuations of the vortices and the quadrupole mode, resulting from the transfer of energy from the vortices to the BEC through excited phonons. When the two quantized vortices are initially positioned sufficiently far from the center, their radial oscillations resonate with the quadrupole mode in the xy plane due to the coincidence of their periods. The kinetic energy components reveal that the quantized vortices experience energetic damping, leading to phonon excitation within the BEC. Consequently, the two rotating quantized vortices spiral out and the quadrupole mode is amplified. This energy transfer does not occur when the two quantized vortices are initially positioned either near the center or near the surface of the BEC. In the former case, this is due to the disparity between the periods of the two oscillations, while in the latter case, it is likely because the corotation of the two quantized vortices is strongly coupled with the quadrupole mode.  \nBased on the findings of this work, we can find that a self-gravitating BEC returns energy to the two corotating quantized vortices after initially receiving energy  \nfrom them. As shown in Fig. 4 (e) and (f), when the initial positions of the two quantized vortices are sufficiently separated, their corotation transitions to an inward spiral after initially spiraling outward to a certain degree. This switching of trajectory indicates a reversal in the direction of energy transfer, with the quantized vortices receiving energy from the BEC. Notably, since the total angular momentum is conserved within this system, we can anticipate that the energy transfer periodically alternates.  \nSuch corotation of two quantized vortices in a selfgravitating BEC can occur following the merger of multiple BECs. When two self-gravitating BECs collide, the resulting BEC often exhibits elliptical deformation accompanied by rotation [54]. Furthermore, as observed in this paper, a cross-sectional view of the elliptical density profile reveals the presence of a few vortex-like structures. These features align with our numerical findings in this study. Therefore, we infer that a self-gravitating BEC after merging is likely to exhibit deformation and rotation driven by collective modes, accompanied by the corotation of two quantized vortices. This post-merger rotation of the self-gravitating BEC and the corotation of the quantized vortices within it could be connected to turbulence within a self-gravitating BEC [37, 38] and may be related to the characteristic rotation of galaxies implying the existence of dark matter [52, 55].  \nThis work leaves some open questions for the corotation of quantized vortices in a self-gravitating BEC. First, while we change only the initial position of the quantized vortices, the dependencies on the strength of the gravitational interaction have yet to be revealed. The strength of the gravitational potential should influence the rotational period of quantized vortices in a self-gravitating BEC by analogy with that in a dipolar BEC. Our result suggests that phonons mediate the energy exchange between the radial oscillations of the quantized vortices and the quadrupole oscillations of the self-gravitating BEC. Thus, similar to the excitation of Kelvin waves in trapped atomic BECs [56], a certain type of phonon decay (or synthesis) may occur.', 'ACKNOWLEDGMENTS': 'This work was supported by JST, the establishment of university fellowships towards the creation of science technology innovation, Grant Number JPMJFS2138. M. T. acknowledges the support from JSPS KAKENHI Grant Number JP23K03305 and JP22H05139.', 'Appendix A: ANALYTICAL CALCULATION OF TWO PARALLEL QUANTIZED VORTICES IN A CIRCULAR CONTAINER BY VORTEX POINT MODEL': "We employ a vortex-point model with friction as a simplified approach to examine the dynamics of two corotating quantized vortices. This model is derived from the vortex filament model, commonly used in the study of quantized vortices in superfluid 4 He [57]. Assuming that the vortex core is sufficiently small, the quantized vortex can be regarded as a filament parametrized by ξ . The vortex filament moves within the superfluid flow field v s ( r , t ), and as a result it obeys the equation of motion  \nd s dt = v s -α s ' × v s + α ' s ' × ( s ' × v s ) , (A1)  \nusing the position of the filament s ( ξ, t ) with the onedimensional coordinate ξ and s ' ( ξ, t ) ≡ ∂ s ( ξ,t ) ∂ξ . The parameters α and α ' imply the strength of dissipation and in the case of superfluid 4 He, they give rise to phenomenological friction between the superfluid and normal-fluid components, so-called 'mutual friction'. However, we assume these parameters do not vanish in Eq. (A1) despite no normal fluid component in order to compare the result to our numerical findings on a self-gravitating BEC. In other words, these dissipations are caused by different mechanisms from the superfluid 4 He. For convenience, we neglect α ' below that.  \nIn an incompressible quantum fluid, the vortex filament L induces the superfluid velocity field v s ( r , t ) following the Biot-Savart law written by  \nv s ( r , t ) = κ 4 π ∫ L dξ s ' ( ξ, t ) × ( r -s ( ξ, t )) | r -s ( ξ, t ) | 3 . (A2)  \nSpecifically, a straight vortex line with circulation κ generates a counter-clockwise velocity field given by  \nv s = κ 2 πr (A3)  \naround the vortex line, where r denotes the distance from the vortex core. Assuming the vortex lines are straight in a cylindrically symmetric system, the system can be reduced to two dimensions, where the vortex lines are treated as points.  \nConsidering quantized vortices are in a container, the boundary effect can be replaced by introducing an image vortex positioned at the corresponding location. In this work, we specifically treat the scenario in which the quantized vortices are confined within a circular container of radius R . The circular container imposes boundary conditions such that the radial component of the velocity field vanishes. Hence, when a vortex point is situated at a distance r v from the center of the container, the image vortex, carrying a charge of -κ , is positioned in the same direction as the quantized vortex, with a distance from the center of the container given by R 2 /r v .  \nFIG. 13. The rotational period of the two quantized vortices corotating in a circular container in Eq. (A5) when R = 6 . 46 kpc and κ = 4 . 10 kpc 2 / Myr. The red dotted line shows the radius where the period becomes maximum and the blue dotted line indicates this value.  \n<!-- image -->  \nWe consider that two vortex points with the same charges are symmetrically located in the container. While one of them is positioned far from the center of the container with a distance r v and an angle θ , the other is positioned with the common distance in the opposite direction. Then, a vortex point moves in the flow field generated by the other and two image vortex points. Using Eqs. (A1) and (A3), the equation of motion is written by  \ndr v dt = ακ 2 π { 1 2 r v + r v R 2 -r 2 v -r v R 2 + r 2 v } , (A4a)  \nr v dθ dt = κ 2 π { 1 2 r v + r v R 2 -r 2 v -r v R 2 + r 2 v } . (A4b)  \nWithout dissipation of the quantized vortices, namely α = 0, we find that the vortex points do not displace radially from the initial positions, as dr v dt = 0, derived from Eq. (A4a). Equation (A4b) then suggests that the vortex points rotate with a constant period  \nT rot = 8 π 2 r 2 o κ R 4 -r 4 o R 4 +3 r 4 o (A5)  \naround the center of the container, where the initial distance of r v is represented as r o . Consequently, the two vortex points undergo uniform circular motion along the same orbit. The r o -dependence of this rotational period is shown in Fig. 13. In this figure, κ = 4 . 10 kpc 2 / Myr is set to the same value as that used for the self-gravitating BEC in this work and R = 6 . 46 kpc corresponds to the Thomas-Fermi radius of this BEC to facilitate comparison to our numerical result. When the vortex points are near the center of the container, the period increases quadratically with r o , similar to that in the uniform case  \ngiven by  \nT rot = 8 π 2 r 2 o κ (A6)  \nwith R → ∞ in Eq. (A5). However, as r o further increases, the period no longer shows sufficient increase and gradually decreases near the boundary due to the influence of the image vortices, i.e. , the finite-size effects. As a reslt, the period reaches a maximum value T rot = √ 2 -√ 3 3 √ 3 8 π 2 κ R 2 ≈ 182 Myr when r o = { 2 √ 3 -1 } 1 4 R ≈ 4 . 05 kpc, as indicated by the red and blue dotted lines.  \nWhen the dissipation of the quantized vortices is present, namely α > 0, the two vortex points move outward while rotating within the container. Solving Eq. (A4a), r v satisfies  \nακ 2 π t = 4 √ 3 9 R 2 tan -1 ( √ 3 r 2 v R 2 ) -4 √ 3 9 R 2 tan -1 ( √ 3 r 2 o R 2 ) -r 2 v -r 2 o 3 , (A7)  \n- [1] A. L. Fetter, Rev. Mod. Phys. 81 , 647 (2009).\n- [2] J. R. Abo-Shaeer, C. Raman, J. M. Vogels, and W. Ketterle, Science 292 , 476 (2001).\n- [3] M. Tsubota, K. Kasamatsu, and M. Ueda, Phys. Rev. A 65 , 023603 (2002).\n- [4] M. C. Tsatsos, P. E. S. Tavares, A. Cidrim, A. R. Fritsch, M. A. Caracanhas, F. E. A. dos Santos, C. F. Barenghi and V. S. Bagnato, Phys. Rep. 622 , 1 (2016).\n- [5] M. Tsubota, M. Kobayashi and H. Takeuchi, Phys. Rep. 522 , 191 (2013).\n- [6] E. P. Gross, J. Math. Phys. Rev. 4 , 195 (1963).\n- [7] L. P. Pitaevskii, Sov. Phys. JETP 13 , 451 (1961).\n- [8] L. Pitaevskii and S. Stringari, Bose-Einstein Condensation , (Oxford University Press, 2003).\n- [9] S. Inouye, M. R. Andrew, J. Stenger, H. -J. Miesner, D. M. Stamper-Kurn and W. Ketterle, Nature 392 , 151 (1998).\n- [10] C. J. Pethick and H. Smith, Bose-Einstein Condensation in Dilute Gases, second edition, (Cambridge University Press, 2008).\n- [11] E. Lundh, A. Collin, and K.-A. Suominen, Phys. Rev. Lett. 92 , 070401 (2004).\n- [12] T. Lahaye, C. Menotti, L. Santos, M. Lewenstein, and T. Pfau, Rep. Prog. Phys. 72 , 126401 (2009).\n- [13] L. Chomaz, I. Ferrier-Barbut, F. Ferlaino, B. LaburtheTolra, B. L. Lev, and T. Pfau, Rep. Prog. Phys. 86 , 026401 (2023).\n- [14] M. Klawunn, R. Nath, P. Pedri, and L. Santos, Phys. Rev. Lett. 24 , 240403 (2008).\n- [15] A. Cidrim, F. E. A. Dos Santos, E. A. L. Henn, and T. Macr'ı, Phys. Rev. A 98 , 023618 (2018).\n- [16] L. Jia, A.-B. Wang, and S. Yi, Phys. Rev. A 97 , 043614 (2018).\n- [17] S. Sabari, R. K. Kumar, and L. Tomio, Phys. Rev. A 109 , 023313 (2024).\n- [18] W. Hu, R. Barkana and A. Gruzinov, Phys. Rev. Lett.  \nwhich demonstrates that r v monotonically increases with time until the vortex points reach the boundary of the container. Furthermore, Eq. (A4b) indicates that the angular velocity also varies over time, initially decreasing and accelerating as the vortex points approach the boundary.  \n85 , 1158 (2000).  \n- [19] T. Matos and L. A. Ure˜na-L'opez, Phys. Rev. D 63 , 063506(2001).\n- [20] C. G. Bohmer and T. Harko, JCAP 06 , 025(2007).\n- [21] P. H. Chavanis, Phys. Rev. D 84 , 043531(2011).\n- [22] E. G. M. Ferreira, Astron. Astrophys. Rev. 29 , 7 (2021).\n- [23] J. C. Niemeyer, Progress in Particle and Nuclear Physics 113 ,103787 (2020).\n- [24] L. Hui, Annu. Rev. Astron. Astrophys. 59 , 247 (2021).\n- [25] T. Matos, L. A. Ure˜na-L'opez, and J.-W. Lee, Front. Astron. Space Sci. 11 , 1347518 (2024).\n- [26] F. S. Guzm'an and L. A. Ure˜na-L'opez, Astrophys. J. 645 , 814 (2006).\n- [27] D. O'Dell, S. Giovanazzi, G. Kurizki, and V. M. Akulin, Phys. Rev. Lett. 84 , 5687 (2000).\n- [28] A. Paredes, D. N. Olivieri, and H. Michiel, Physica D: Nonlinear Phenomena 403 , 132301 (2020).\n- [29] S. Giovanazzi, G. Kurizki, I. E. Mazets, and S. Stringari, Europhys. Lett., 56 , 1 (2001).\n- [30] K. Asakawa, H. Ishihara, and M. Tsubota, Prog. Theor. Exp. Phys. 6 , 063J01 (2024).\n- [31] Y. O. Nikolaieva, Y. M. Bidasyuk, K. Korshynska, E. V. Gorbar, J. Jia, and A. I. Yakimenko, Phys. Rev. D 108 , 023503 (2023).\n- [32] N. T. Zinner, Phys. Res. Int. 2011 , 734543 (2011).\n- [33] B. Kain and H. Y. Ling, Phys. Rev. D 82 , 064042 (2010).\n- [34] T. Rindler-Daller and P. R. Shapiro, MNRAS 422 , 135 (2012).\n- [35] Y. O. Nikolaieva, A. O. Olashyn, Y. I. Kuriatnikov, S. I. Vilchinskii, and A. I. Yakimenko, Low Temperature Physics 47 , 684 (2021).\n- [36] A. S. Dmitriev, D. G. Levkov, A. G. Panin, E. K. Pushnaya, and I. I. Tkachev, Phys. Rev. D 104 , 023504 (2021).\n- [37] P. Mocz, M. Vogelsberger, V. H. Robles, J. Zavala, M. Boylan-Kolchin, A. Fialkov, and L. Hernquist, MNRAS\n- 471 , 4559 (2018).\n- [38] I-K. Liu, N. P. Proukakis, and G. Rigopoulos, MNRAS 521 , 3625 (2023).\n- [39] C. F. Barenghi, H. A. J. Middleton-Spencer, L. Galantucci, and N. G. Prker, AVS Quantum Science 5 , 025601 (2023).\n- [40] R. Navarro, R. Carretero-Gonz'alez, P. J. Torres, P. G. Kevrekidis, D. J. Frantzeskakis, M. W. Ray, E. Altunta¸s, and D. S. Hall, Phys. Rev. Lett. 110 , 225301 (2013).\n- [41] S. Gautam, A. Roy, and S. Mukerjee, Phys. Rev. A 89 , 013612 (2014).\n- [42] K. Nakamura, D. Babajanov, D. Matrasunov, M. Kobayashi, and P. Muruganandam, J. Phys. A: Math. Theor. 49 , 315102 (2016).\n- [43] Q. Zhao, J. Low Temp. Phys. 204 , 1 (2021).\n- [44] Q. Zhao, H.-J. Bi, X.-M. Yang, and L.-L. Zhang, Acta Phys. Pol. A 141 , 578 (2022).\n- [45] C. Nore, M. Abid and M. E. Brachet: Phys. Rev. Lett. 78 , 3896 (1997).\n- [46] C. Nore, M. Abid and M. E. Brachet: Phys. Fluids 9 , 2644 (1997).\n- [47] F. S. Guzm'an, F. D. Lora-Clavijo, J. J. Gonz'alez-Avil'es and F. J. Rivera-Paleo, Phys. Rev. D 89 , 063507 (2014).\n- [48] Y. Nesterov, Soviet Math. Dokl. 27, 372 (1983).\n- [49] B. O'Donoghue and E. J. Cand'es, Found. Comput. Math. 15 , 715 (2015).\n- [50] W. Su, S. Boyd and E. J. Cand'es, J. Mach. Learn. Res. 17 , 153 (2016).  \n- [51] M. Abad, M. Guilleumas, R. Mayol, M. Pi, and D. M. Jezek, Phys. Rev. A 79 , 063622 (2009).\n- [52] J. Binney and S. Tremaine, GALACTIC DYNAMICS, second edition, (Princeton University Press, 2008).\n- [53] T. K. Ghosh, Phys. Rev. A 65 , 053616 (2002).\n- [54] B. Schwabe, J. C. Niemeyer, and J. F. Engels, Phys. Rev. D 94 , 043513 (2016).\n- [55] V. Rubin and W. Ford, Astrophys. J, 159 , 379 (1970).\n- [56] T. Mizushima, M. Ichioka, and K. Machida, Phys. Rev. Lett. 90 , 180401 (2003).\n- [57] R.J. Donnelly, Quantized Vortices in Helium II , (Cambridge University Press, Cambridge, England, 1991).\n- [58] See Supplemental Material for the time evolution of the mass density profile in the plane z = 0 when r o ≈ 2 . 0 kpc.\n- [59] See Supplemental Material for the time evolution of the vortex cores in the three-dimensional system when r o ≈ 1 . 0 kpc.\n- [60] See Supplemental Material for the time evolution of the mass density profile in the plane z = 0 when r o ≈ 3 . 0 kpc.\n- [61] The reduction in the size of the BEC is attributed to stronger gravitational binding resulting from the increased central mass density. When two quantized vortices are positioned far from the center of the BEC, they do not interfere with the accumulation of particles in the central region. Consequently, the central density rises with increasing r o , thereby enhancing the gravitational potential."}