Self-interacting multistate boson stars

In this paper, we consider rotating multistate boson stars with quartic self-interactions. In contrast to the nodeless quartic-boson stars in \cite{Herdeiro:2015tia}, the self-interacting multistate boson stars (SIMBSs) have two types of nodes, including the $^1S^2S$ and $^1S^2P$ states. We show the mass $M$ of SIMBSs as a function of the synchronized frequency $\omega$, and the nonsynchronized frequency $\omega_2$ for three different cases. Moreover, for the case of two coexisting states with self-interacting potential, we study the mass $M$ of SIMBSs versus the angular momentum $J$ for the synchronized frequency $\omega$ and the nonsynchronized frequency $\omega_2$. Furthermore, for three different cases, we analyze the coexisting phase with both the ground and first excited states for SIMBSs. We also calculate the maximum value of coupling parameter $\Lambda$, and find the coupling parameter $\Lambda$ exists the finite range.


I. INTRODUCTION
The scalar field dark matter (SFDM) is one of the important models of the dark matter [2][3][4][5], and the model considers a spin-0 scalar field of the very small mass (typically ∼ 10 −22 eV/c 2 ). Due to the ultra-light boson mass, the SFDM could form Bose-Einstein condensates (BEC) in the very early Universe, which are interpreted as the dark matter haloes. The presence of dark matter halo is inferred from its gravitational effect on a spiral galaxy rotation curve (RC). At this point boson stars (BSs) could play an important role.
The BSs were studied in the pioneering works of Kaup et al. [6,7], where they considered a complex scalar field coupled to Einstein gravity theory, ones refer to this coupled system as the Einstein-Klein-Gordon (EKG) theory. See Refs. [8,9] for a review.
For the spherical symmetry solutions of boson stars, Deng and Huang [10] have constructed a system consisting of coexisting states of two scalar fields, including two ground states, and they calculated the gravitational redshift of the system. In 2010, Bernal et al. [11] have obtained the configurations with two states, a ground and first existed states, as well as they have produced the RC of multistate boson stars that are flatter at large radii than the RC of boson stars with only one state. The study of spherical symmetry solutions of boson stars can be extended to the gauge charge [12][13][14] case, with de Sitter [15,16], or anti-de Sitter [17,18] boundary conditions and the geodesics case [19]. Moreover, the axially symmetric solutions of boson stars were studied by Schunck and Mielke [20,21], Yoshida and Eriguchi [22], Herdeiro and Radu [23]. However, they considered only ground state of the free, complex scalar field. The study of rotating axisymmetric solutions of BSs can be extended to the excited case [24], the multistate boson stars [25].
In 1986, Colpi et al. [26] have constructed the self-interacting boson stars solutions with the type λ|ψ| 4 , the sixth power |ψ|-terms potentials were obtained by Mielke and Scherzer [27] and Kleihaus et al. [28]. Generalizing the rotating Q-balls models of Volkov and Wöhner [29] to the self-gravitating case. In addition, in 2015, Herdeiro and Radu have studied the quartic-BSs solutions [1], which is extended to the Kerr black holes with selfinteracting scalar hair. The study of BSs can be extended to the case of charged BSs [30][31][32], nonrotating and rotating axion boson stars [33,34], Proca stars [35][36][37], and in the higherdimensional spacetime case [38]. The linear stability of boson stars with respect to small oscillations was discussed by Lee and Pang in [39], the study of the stability properties of boson stars extended to the quartic and sextic self-interaction term case [40], dynamical -BSs [41,42], and non-linear dynamics of spinning BSs case [43]. The properties of BSs have been investigated widely [44][45][46][47][48][49][50][51][52][53][54]. In the present work, we would like to numerically solve the EKG equations and give a family of rotating multistate boson stars for three different cases, including the ground state with self-interacting potential, the first excited state with self-interacting potential and the two coexisting states with self-interacting potential.
The paper is organized as follows. In Sect. II, we introduce the model of the fourdimensional Einstein gravity coupled to two self-interacting complex scalar fields ψ i (i = 1, 2). In Sect. III, the boundary conditions of the SIMBSs are studied. We exhibit the numerical results of the equations of motion and show the properties of the 1 S 2 S and the 1 S 2 P states for three different cases in Sect. IV. Conclusions and perspective are given in Sect. V.

II. THE MODEL SETUP
We consider two massive, complex scalar fields model coupled minimally to gravity in 3 + 1 dimensions describing self-interacting multistate boson stars. The action is where R is the Ricci scalar, G is Newton's constant and L m denotes the matter Lagrangian: where λ, µ i , (i = 1, 2) are the interaction parameter and the scalar fields mass, respectively.
It is convenient to change the dimensionless coupling parameter λ to where M 2 P l = G −1 is the Planck mass. For simplicity of presentation, we will set G = c = µ 1 = 1. The corresponding equations of motion are given by where is the covariant D'Alembert differential operator.
To obtain stationary axisymmetric solutions of self-interacting multistate boson stars, we choose the ansatz as follows, see e.g. [22,23,25]: In addition, the ansatz of two complex scalar fields ψ i are given by Here F 0 , F 1 , F 2 , W , and φ i(n) (i = 1, 2) are functions of the radial distance r and the polar angle θ, only. Subscript n of Eq. (7) is named as the principal quantum number of the scalar field, and n = 0 is regarded as the ground state and n ≥ 1 as the excited states. Subscript i are indicated by two complex scalar fields only. In additions, the constants m i (i = 1, 2) and ω i (i = 1, 2) are indicated by the azimuthal harmonic index and the frequency of the complex scalar field, respectively. When ω 1 = ω 2 = ω, the frequency of the scalar field is called the synchronized frequency, while ω 1 = ω 2 is called the nonsynchronized frequency.
It is well known the rotating boson stars with first excited state exists two types of nodes, including radial and angular nodes. The radial and angular nodes can seen in middle panels of Figs. 1, 4 of [25].
At the origin we require We note that the values of F 0 (0, θ), F 1 (0, θ), F 2 (0, θ), and W (0, θ) are the constants that are not dependent of the polar angle θ.
Near the boundary r → ∞, on the other hand, the mass of boson stars M and total angular momentum J are extracted from the asymptotic behavior of the metric functions,

IV. NUMERICAL RESULTS
In this section, we will solve the above coupled Eqs. (5a) and (5b) with the ansatzs (6) and (7) numerically. It is convenient to introduce a new coordinate x ≡ r 1+r , which implies that the new radial coordinate x ∈ [0, 1]. Thus, the inner and outer boundaries of the shell are fixed at x = 0 and x = 1, respectively. By exploiting the reflection symmetry θ → π − θ on the equatorial plane, it is enough to consider the range θ ∈ [0, π/2] for the angular variable. All numerical calculations are based on the finite element methods. Typical grids used have sizes of 100×100 in the integration region 0 ≤ x ≤ 1 and 0 ≤ θ ≤ π 2 . Our iterative process is the Newton-Raphson method, and the relative error for the numerical solutions in this work is estimated to be below 10 −5 .
Next, we will discuss the SIMBSs, including the principal quantum number n = 0, which is the ground state case and the principal quantum number n = 1, which belongs to the case of the first excited state. Besides, we exhibit two classes of radial n r = 1 and angular n θ = 1 node solutions, respectively. As noted above, for the case of rotating boson stars with first excited state, there is a similar situation in atomic theory and quantum mechanics, the first excited state of hydrogen has an electron in the 2s orbital and 2p orbital, which correspond to the radial and angle node, respectively. Therefore, the coexisting states of two scalar fields, which have a ground state and a first excited state with a radial node n r = 1, is named as the 1 S 2 S state as well as the coexistence of a ground state and a first excited state with an angle node n θ = 1 is called the 1 S 2 P state.

1 S 2 S state
In this subsection, we will analyze the solutions with an even-parity scalar field. The value of the scalar field φ i (i = 1, 2) can be seen in of Fig. 1 of [25]. According to the numerical results, we find that there exists a maximal mass of multistate boson stars for parameter Λ 1 < Λ c = 2.8. To study the properties of the SIMBSs, we exhibit in Fig In the left panel of Fig. 1, we exhibit the mass M versus the synchronized frequency ω with the azimuthal harmonic index m 2 = 1 , 2, 3, denoted by the origin, blue and red lines, respectively, and the shaded regions indicate 0 ≤ Λ 1 ≤ 2. More details are given in the inset plotted in the left panel of Fig. 1. We observe that for the same synchronized frequency and coupling parameters, the free, multistate boson stars possess a higher mass than the SIMBSs. Moreover, we observe that for the fixed coupling parameter Λ 1 = 0 and Λ 1 = 2, the multistate boson stars exist only a stable branch, which is different from the quartic-BSs case in Ref. [1]. In the right panel of Fig. 1, we show the mass M versus the nonsynchronized frequency ω 2 for the fixed value of ω 1 = 0.85, as well as the shaded regions indicate 0 ≤ Λ 1 ≤ 2. We can see that the free, multistate boson stars presents a higher value for the mass, which is different from the quartic-BSs in Ref. [1]. From the inset plotted of the right of Fig. 1, we found that the minimum value of mass of the SIMBSs is  increases, the domain of existence of µ 2 decreases. Here we note that, the SIMBSs case does not occur with another branch. Moreover, we find the another branch do not depend on the values of λ, in order to verify whether there exists another family of multistate solutions, we also use the same method in [25] to seek for the new family of SIMBSs case (The more detailed analysis process can be seen from the sixth paragraph of section IV (A) and in Fig.   3 of [25]).

1 S 2 P state
In order to compare the results with even-parity scalar field, in this subsection we also study the solutions with odd-parity scalar field mentioned in [1,24,25]. Meanwhile, in the middle panels of Fig. 5 of [25], we can see that the odd function φ 2 is antisymmetric with respect to the x-axis for the angular variable 0 ≤ θ ≤ π 2 . Besides, according to the numerical results, we find no solutions of self-interacting multistate boson stars for Λ 1 ≥ Λ c = 5.8. (origin shaded area), m 2 = 2 (blue shaded area), m 2 = 3 (red shaded area), respectively. In the inset panel, we show the detail of the blue shaded area with m 2 = 2. We observe that the SIMBSs exhibits only a stable branch with m 2 = 1 , 2, 3, which is similar to the case of free, multistate boson stars in Ref. [25]. We again observe that with the decrease of the frequency ω, the mass of SIMBSs keeps increasing. However, we can see that the 1 S 2 P state curves exhibit only a stable branch with m 2 = 1 , 2, 3, which is similar to the case of free, multistate boson stars in Ref. [25], as well as the mass of the SIMBSs is smaller than the case of free, multistate boson stars.
In the right panel of Fig. 2, we plot the mass of the SIMBSs versus the nonsynchronized frequency ω 2 for the ω 1 = 0.85. We observe that the free, multistate boson stars have higher mass level than the SIMBSs. In Table II, we present the domain of existence of the mass µ 2 of the scalar field φ 2 for the different values of Λ 1 . We note that with the increase of coupling parameter Λ 1 , the domain of existence of the µ 2 begins to decrease.

1 S 2 S state
According to the numerical results, we find that the range of coupling parameter Λ 2 ∈ (0, 388). In the left panel of Fig. 3, we show the mass M of the SIMBSs versus the synchronized frequency ω with the azimuthal harmonic index m 2 = 1 , 2, 3, represented by the origin, blue, and red lines, respectively, and the shaded regions indicate 0 ≤ Λ 2 ≤ 10. More details are given in the inset plotted in the left panel of Fig. 3. We observe that for the fixed coupling parameters Λ 2 = 0 and Λ 2 = 10, all the multistate boson stars curves exist a unique branch, which is different from the case of ground state in Ref. [1]. Moreover, we observe that for the same synchronized frequency and coupling parameters, the self-interacting multistate boson stars solutions is heavier than the case of free, multistate boson stars. This result can be understood that when the synchronized frequency ω tends to its maximum, the scalar field of the first excited state could reduce to zero and there exists only a single scalar field of the ground state. So, three sets of the multistate boson stars intersect with the ground state solutions, respectively. In the right panel of Fig. 3, we present the mass M of the SIMBSs versus the nonsynchronized frequency ω 2 for the fixed value of ω 1 = 0.85, and the dashed lines and solid lines stand for Λ 2 = 0 and Λ 2 = 10, respectively. From the graphic, it is obvious that the SIMBSs has higher mass, which is similar to the case of quartic-BSs in Ref. [1]. In the inset plotted of the right of

1 S 2 P state
As typical examples for our numerical results, we obtain the range of coupling parameter Λ 2 ∈ (0, 41). Meanwhile, we show in Fig. 4 the SIMBSs mass M as functions of the synchronized frequency ω and nonsynchronized frequency ω 2 , respectively, by taking m 2 = 1 (origin shaded region), m 2 = 2 (blue shaded region), m 2 = 3 (red shaded region). In the left panel, we observe that as azimuthal harmonic index m 2 increases, the maximum value of the synchronized frequency ω increases, and the minimum value of the mass increases as well.
In the right panel of Fig. 4 In the above two subsections, we gave two families of self-interacting multistate boson stars. Next, we consider the self-interacting potential U (|ψ 1 |) = 0, U (|ψ 2 |) = 0. The matter m is given by with U (|ψ i |) = µ 2 i |ψ i | 2 + λ i |ψ i | 4 , we still use the same relation between the parameters Λ i and λ i given in Eq.(4), i = 1, 2.
We can see that the multistate boson stars curves have the similar behavior as the above two subsections case of 1 S 2 S state. With the increase of the synchronized frequency ω, the mass M of the multistate boson stars decreases. In the right panel of Fig. 5, we show the mass M of the multistate boson stars vs the nonsynchronized frequency ω 2 with fixed value of ω 1 = 0.85 for three sets of the azimuthal harmonic index. We found that the mass M of the self-interacting multistate boson stars is larger than that of the free, multistate boson stars, under the common influence, which are both the coupling constants Λ 1 and Λ 2 . Besides, in the inset of the right of Comparing with the results of the ground state boson stars in Ref. [1], we can see that the case of the SIMBSs does not occur in zigzag patterns, and the minimum value of the mass M is larger than that of the ground state boson stars. In order to explore the influence  Table V. We observe that the range of mass µ 2 decreases with the increase of the coupling parameters Λ 1 and Λ 2 , which is similar to the case of 1 S 2 S state in above two subsections.

1 S 2 P state
According to the numerical results, we obtain the range of coupling parameters Λ 1 ∈ (0, 5.8), Λ 2 ∈ (0, 41). In the left panel of Fig. 7, we plot the mass M of the SIMBSs vs the synchronized frequency ω with the azimuthal harmonic index m 2 = 1 (the origin area), m 2 = 2 (the blue area), m 2 = 3 (the red area), respectively. In the inset of the left panel of Fig. 7, we show the detail of the blue and red shaded regions with m 2 = 2 and m 2 = 3.
From the graphic, we found that the SIMBSs exists only a stable branch, which is similar to the 1 S 2 P state in above two subsections. Moreover, we found that the mass of multistate boson stars decreases with increasing the frequency ω. In the right panel of Fig. 7 which is similar to the case of 1 S 2 S state of the two coexisting states with self-interacting potential case.

V. CONCLUSION
In this paper, we have constructed rotating boson stars composed of the coexisting states of two quartic-self-interacting scalar fields, including the ground and first excited states.
Comparing with the solutions of the quartic-BSs in Ref. [1], we have found that the SIMBSs have two types of nodes, including the 1 S 2 S and 1 S 2 P states. According to the numerical results, we find that the coupling parameter Λ is the finite value for the 1 S 2 S and 1 S 2 P states, which is different from the quartic-BSs in Ref. [1]. By calculating the coexisting phase of the SIMBSs, we found that the domain of existence of the mass µ 2 decreases with increasing the coupling parameters. When the synchronized frequency ω tends to its maximum, the first excited state could reduce to zero, and there exists only the ground state. However, for the case of nonsynchronized frequency ω 2 , we fixed the value of ω 1 = 0.85, the minimal mass of the SIMBSs is always a constant value for the different azimuthal harmonic indexes m 2 = 2, 3.
It is important to understand how the stability properties of the self-interacting multistate boson stars, according to numerical analysis of the stability of the self-interaction boson stars in [40], the authors found that the solutions have both stable and unstable branches. Besides, in 2010, Bernal et al. [11] found for the nonrotating multistate boson stars solutions, such configurations are stable if the number of particles in the ground state is larger than the number of particles in the excited state. To our knowledge, it is difficult to numerically analyze the stability properties of self-interaction multistate boson stars. However, a good way to guarantee the stability is to have the scalar field perturbations decouple from the metric perturbations in Ref. [56].
There are several interesting extensions of our work. Firstly, we have studied the SIMBSs, next, we will investigate the gravitational redshift and the RC of multistate boson stars inspired by the work [10,57]. In 2019, Brihaye and Hartmann have studied the spontaneous scalarization of spherically symmetric, asymptotically flat boson stars in Ref. [58].
Therefore, the second extension of our work is to study the spontaneous scalarization for axially symmetric BSs solutions. Finally, we are planning to construct this solutions where the self-interaction potential of the scalar field ψ coupled minimally to Einstein gravity in higher-dimensional spacetime.