A Newtonian Analysis of Multi-scalar Boson Stars with Large Self-couplings

TL;DR

This paper analyzes multi-scalar boson stars with large self-couplings using Newtonian approximation, aiming to model exotic mass distributions that could explain galaxy rotation curves.

Contribution

It introduces a novel approach linking coupling matrices to graph theory to construct models with unique mass distributions for boson stars.

Findings

Models produce exotic mass distributions.

Coupling matrices related to graph structures.

Potential explanations for galaxy rotation curves.

Abstract

We study solutions for boson stars in the multi-scalar field theory with global symmetry $[U(1)]^N$. The properties of the boson stars are investigated by the Newtonian approximation with the large coupling limit. Our purpose is to study the models bringing about exotic mass distributions which explain flat rotation curves of galaxies. We propose plausible models in which coupling matrices are associated with various graphs in graph theory.

Equations

$$S$$

$$R_{\nu }^{\mu }-\frac{1}{2}{\delta }_{\nu }^{\mu }R=8\pi G{T}_{\nu }^{\mu },$$

$$T_{\mu\nu}=\sum_{i}\Big{[}\partial_{\mu}\phi_{i}^{*}\partial_{\nu}\phi_{i}+\partial_{\nu}\phi_{i}^{*}\partial_{\mu}\phi_{i}-g_{\mu\nu}|\partial\phi_{i}|^{2}-g_{\mu\nu}m^{2}|\phi_{i}|^{2}\Big{]}-g_{\mu\nu}V(|\phi_{i}|^{2})\,,$$

$$V(|{\varphi }_i{|}^2)\equiv \frac{1}{2}\sum _{i,j}{\tilde{\lambda }}_{ij}|{\varphi }_i{|}^2|{\varphi }_j{|}^2.$$

$$\square {\varphi }_i-m^2{\varphi }_i-\sum _j{\tilde{\lambda }}_{ij}|{\varphi }_j{|}^2{\varphi }_i=0,$$

$$g_{00}\approx -\left(1+2\Phi \right),g_{ij}\approx {\delta }_{ij},g_{0i}=g_{i0}\approx 0,\sqrt{-g}\approx 1,$$

$$\phi_{i}\approx\frac{1}{\sqrt{2m}}\psi_{i}(\mbox{\boldmath$r$},t)\,e^{-imt}\,,$$

$$i{\dot{\psi }}_i=-\frac{1}{2m}{\nabla }^2{\psi }_i+m\Phi {\psi }_i+\frac{1}{4}\sum _j\frac{{\tilde{\lambda }}_{ij}}{m^2}|{\psi }_j{|}^2{\psi }_i$$

$${\nabla }^2\Phi =4\pi Gm\sum _i|{\psi }_i{|}^2.$$

$$S\cong\int dt\,d^{3}\mbox{\boldmath$r$}\left[-\frac{1}{8\pi G}({\nabla}\Phi)^{2}+\sum_{i}\left\{i\psi_{i}^{*}\dot{\psi}_{i}-\frac{1}{2m}|\nabla\psi_{i}|^{2}-m\Phi|\psi_{i}|^{2}\right\}-\frac{1}{8}\sum_{i,j}\frac{\tilde{\lambda}_{ij}}{m^{2}}|\psi_{i}|^{2}|\psi_{j}|^{2}\right]\,,$$

$$\displaystyle\hat{H}(\psi_{i},\Phi)=\int d^{3}\mbox{\boldmath$r$}\,{\cal H}$$

$$\hat{n}_{i}\equiv\int d^{3}\mbox{\boldmath$r$}\,|\psi_{i}|^{2}\,.$$

$${\nabla }^2\Phi =4\pi Gm\sum _i|{\psi }_i{|}^2,$$

$$-\frac{1}{2m}{\nabla }^2{\psi }_i+m\Phi {\psi }_i+\frac{1}{4}\sum _j\frac{{\tilde{\lambda }}_{ij}}{m^2}|{\psi }_j{|}^2{\psi }_i={\mu }_i{\psi }_i.$$

$$\Phi \to \Phi +{\Phi }_0,{\mu }_i\to {\mu }_i+m{\Phi }_0.$$

$${\Lambda }_{ij}\equiv \frac{1}{8\pi G{m}^2}{\tilde{\lambda }}_{ij}.$$

$$\mbox{\boldmath$r$}_{*}\equiv\frac{m}{\sqrt{\Lambda}}\mbox{\boldmath$r$}\,,\quad\Psi_{i}\equiv\sqrt{\frac{4\pi G\Lambda}{m}}\psi_{i}\,,\quad{u_{i}}\equiv\frac{\mu_{i}}{m}\,,$$

$${\nabla }_{\ast }^2\Phi =\sum _i|{\Psi }_i{|}^2,$$

$$-\frac{1}{2\Lambda }{\nabla }_{\ast }^2{\Psi }_i+\Phi {\Psi }_i+\frac{1}{2}\sum _j{C}_{ij}|{\Psi }_j{|}^2{\Psi }_i=u_i{\Psi }_i,$$

$$\Phi {\Psi }_i+\frac{1}{2}\sum _j{C}_{ij}|{\Psi }_j{|}^2{\Psi }_i=u_i{\Psi }_i.$$

$$\Phi \to {\lambda }^2\Phi ,u_i\to {\lambda }^2{u}_i,{\Psi }_i\to \lambda {\Psi }_i,$$

$${\rho }_i\equiv |{\Psi }_i{|}^2.$$

$$n_{i}=\frac{\sqrt{\Lambda}}{4\pi Gm^{2}}\int d^{3}\mbox{\boldmath$r$}_{*}\,\rho_{i}(\mbox{\boldmath$r$}_{*})\,.$$

$$E=\frac{\sqrt{\Lambda}}{4\pi Gm}\int d^{3}\mbox{\boldmath$r$}_{*}\left[\frac{1}{2}({\nabla}_{*}\Phi)^{2}+\Phi\sum_{i}|\Psi_{i}|^{2}+\frac{1}{4}\sum_{ij}C_{ij}|\Psi_{i}|^{2}|\Psi_{j}|^{2}\right]\,.$$

$$E$$

$$M=m\sum _i{n}_i+E=\sum _i{n}_i\left(m+\frac{1}{2}{\mu }_i\right)=m\sum _i{n}_i\left(1+\frac{1}{2}{u}_i\right).$$

$$\frac{1}{x}\frac{d^2}{d{x}^2}(x\Phi (x))=S(x),$$

$$\left(\Phi (x)+\frac{1}{2}\sum _j{C}_{ij}{\rho }_j(x)-u_i\right){\rho }_i(x)=0.$$

$$S(x)\equiv\sum_{i}\rho_{i}(x)\,,\quad x\equiv\sqrt{\mbox{\boldmath$r$}_{*}\cdot\mbox{\boldmath$r$}_{*}}\,.$$

$$\frac{1}{x}\frac{d^2}{d{x}^2}(x\Phi (x))={\rho }_1(x)+{\rho }_2(x),$$

$$\left[\Phi (x)+\frac{1}{2}\left(C_{11}{\rho }_1(x)+C_{12}{\rho }_2(x)\right)-u_1\right]{\rho }_1(x)=0,$$

$$\left[\Phi (x)+\frac{1}{2}\left(C_{21}{\rho }_1(x)+C_{22}{\rho }_2(x)\right)-u_2\right]{\rho }_2(x)=0,$$

Full text

I Introduction

From the perspective of particle cosmology, the boson star Jetzer ; LM ; SM ; LP is one of the candidates for dark matter R1SL1 ; R1SL2 ; LK ; TCL ; ST ; ABBR ; MA ; BBAP ; FG ; Sin ; JS . The ($1+3$)-dimensional boson star was studied as the simplest object of a self-gravitating system and the Newtonian treatment of gravitating bosons has been often discussed (for example, see Refs. Jetzer ; SM ; LP ; LK ; MA ; FG ; Sin ; JS ; KSrecent ).

Boson stars or boson halos have been studied sometimes in expectation of solving the flat rotation curve problem LK ; TCL ; ST ; ABBR ; MA ; BBAP ; FG ; Sin ; JS .111The flat rotation curve problem was also approched by Schunck Schunck , who considered a configuration of a massless scalar field with infinite range. Some authors have attempted to explain rotation curves of galaxies by assuming the existence of a galactic-scale boson star located at the center of a galaxy. In order to fit the observed data, the mass density of the boson star needs to be widely distributed. Such configurations can be constructed by the models such as boson stars with a scalar field in an excited state or rotating boson stars and show better agreement with the galaxy rotation curves, but these boson stars are sometimes unstable. Whereas Newtonian boson stars in only one ground state are stable, it is difficult to illustrate the realistic rotation curve owing to their compact density distribution. Alternative models to solve these problems have been investigated by Matos and Ureña-López MA , and by Bernal et al. BBAP . They considered the multi-state boson stars, i.e. a scalar field both in ground and in excited states, with no (quartic) self-couplings.

In the present paper, we consider scalar boson stars made of multi-species or charges, not multi-states of a single scalar field. We suggest the models which contain multiple scalar fields with self- and mutual-couplings.222Interacting boson stars and Q-balls have been studied by Brihaye et al. Brihaye1 ; Brihaye2 ; Brihaye3 in the other context. It is shown that the configuration induced from the multi-bosons can improve the flatness of the galactic rotation curve at the large scale.

We also claim that the coupling matrix can be built up by the knowledge of graph theory. The model we propose contains $N$ scalar fields interacting with oneself and with ‘adjacent’ scalar fields on a graph. A similar type of many $U(1)$ interacting fields has been motivated by the graph-oriented model with supersymmetry KKS .

To discuss the qualitative behavior of models such as a static and spherical boson star, we study the system of scalar fields with large self-couplings colpi in this paper. An understanding of the basic aspects of general multi-boson systems is yet expected from this approximation.

The present paper is organized as follows. In Sec. II, we propose the action, the Hamiltonian, and the field equations for a model of self-interacting, gravitating scalar fields in the Newtonian limit. The large coupling limit of the model is defined in Sec. III. In Sec. IV, as the simplest case, analytic solutions for boson stars in a bi-scalar theory are obtained, and relations among their physical quantities is derived. In Sec. V, we study boson stars in a multi-scalar theory with only self-couplings of individual scalar fields, i.e., with a diagonal coupling matrix. We propose symmetric models based on the graph theory, which are compatible to deriving the desired profile of a multi-scalar boson star, in Sec. VI. In this section, we mainly concentrate ourselves on a model associated with a complete graph. Other possible models based on the other types of graphs are also mentioned. The last section VII is devoted to a summary and future prospects.

II The multi-scalar model and its Newtonian limit

We consider a system of interacting, gravitating complex scalar bosons $\phi_{i}$ $(i=1,\dots,N)$ of equal mass $m$,333In the large coupling limit being treated in the next section, a possible mass spectrum is absorbed into redefinition of scalar fields by construction. governed by the following relativistic action (where $\hbar=c=1$);

[TABLE]

where d^{4}x=dt\,d^{3}\mbox{\boldmathr} , $G$ is the Newton constant, $R$ is the Ricci scalar, $g$ is the determinant of the metric $g_{\mu\nu}$ $(\mu,\nu=0,1,2,3)$, $\tilde{\lambda}_{ij}=\tilde{\lambda}_{ji}$ $(i,j=1,\dots,N)$ are the dimensionless real (self and mutual) coupling constants, and $|\partial\phi_{i}|^{2}\equiv g^{\mu\nu}(\partial_{\mu}\phi_{i})^{*}(\partial_{\nu}\phi_{i})$. The action has global $[U(1)]^{N}$ symmetry.

By the variational principle, we derive the Einstein equation from the action as

[TABLE]

where the energy-momentum tensor $T_{\mu\nu}$ in the system is given by

[TABLE]

where the quartic interaction is denoted as

[TABLE]

On the other hand, the equation of motion for each scalar field $\phi_{i}$ is given by

[TABLE]

where $\Box\phi\equiv\frac{1}{\sqrt{-g}}\partial_{\mu}(\sqrt{-g}g^{\mu\nu}\partial_{\nu}\phi)$ is the covariant d’Alembertian.

Throughout the present paper, we consider the gravitating system in the Newtonian approximation. The Newtonian limit can be attained by assuming that the spacetime metric in the weak field approximation can be written as

[TABLE]

where $\Phi$ is the Newtonian gravitational potential.

Assuming further that the complex scalar field has a nearly harmonic time dependence expressed by

[TABLE]

we obtain the (non-linear) Schrödinger equation

[TABLE]

as the Newtonian limit of Eq. (5), where $\nabla^{2}$ is the Laplacian in the flat space and the dot ($\dot{~{}}$) indicates the time derivative. In the present limit, the Einstein equations reduce to the Poisson equation

[TABLE]

The Newtonian treatment of the Lagrangian and Hamiltonian is as follows. We find the following Newtonian action in the limit:

[TABLE]

where $(\nabla\Phi)^{2}\equiv{\bf\nabla}{\Phi}\cdot{\bf\nabla}{\Phi}$ and the symbol $\cong$ indicates that some surface terms have been omitted. Therefore, the Hamiltonian of the system is derived as

[TABLE]

The number of particles of the $i$-th species is expressed as

[TABLE]

In addition, we require the condition $\hat{n}_{i}=n_{i}$, i.e., the condition that the system contains $n_{i}$ scalar bosons of the $i$-th species. Then, we consider ${\delta}\{\hat{H}-\sum_{i}\mu_{i}(\hat{n}_{i}-n_{i})\}=0$ as equations for the scalar fields in the mean field approximation, where $\mu_{i}$ are Lagrange multipliers and can be interpreted as chemical potentials for corresponding bosons.

Now, we obtain coupled equations for the stationary gravitational field and the scalar fields as follows:

[TABLE]

Therefore, the system is reduced in the Newtonian limit to the (non-linear) Schrödinger–Poisson system.

It is notable that the field equations (13) and (14) are invariant under a common shift of potentials

[TABLE]

Therefore, we can choose $\Phi=0$ at the spatial infinity for a compact boson star, even after solving the field equations.

In the subsequent sections of this paper, we will restrict ourselves on the large coupling limit for compact objects, which will be defined in the next section.

III The large coupling limit and the field equations

Here, we consider the large coupling limit colpi . It is incidentally known that the large coupling leads to a large scale boson star. First we define the matrix of couplings $\Lambda_{ij}$ as

[TABLE]

In addition, we introduce the following quantities

[TABLE]

where $\Lambda$ is a typical scale of $\Lambda_{ij}$. Then, the set of equations reduces to the simple form

[TABLE]

where $C_{ij}\equiv\frac{1}{\Lambda}\Lambda_{ij}$. $\nabla^{2}_{*}$ is the rescaled Laplacian expressed in terms of the coordinate \mbox{\boldmathr}_{*}.

In the limit of $\Lambda\rightarrow\infty$, equation (19) further reduces to the algebraic equation:

[TABLE]

In this paper, we only consider the case with the large coupling limit. It is interesting to note that the equation (20) as well as the Poisson equation (18) are invariant under the following scale transformation:

[TABLE]

where $\lambda$ is a constant.

We further define normalized (particle number) density functions as

[TABLE]

Then, the particle number of the $i$-th scalar boson is given by

[TABLE]

The Newtonian energy $E$ of the system in the large coupling limit can be expressed from Eq. (11) as

[TABLE]

Substituting the solution of Eqs. (18), (20) and (22) into this equation, we obtain

[TABLE]

For a compact object, the energy becomes negative ($E<0$).

The mass of the boson star is given in the present Newtonian scheme by

[TABLE]

Finally in the present section, we consider the field equations for static, spherically symmetric solutions of the system. Then, the field equations (18) and (19) can be rewritten as

[TABLE]

[TABLE]

where

[TABLE]

IV Non-Relativistic spherical bi-scalar boson stars

in large coupling limit

We first examine the simplest case, static and spherical solutions for boson stars in the system with two scalar fields. There are not single, but two scalar conserved charges for boson 1 and boson 2.

We consider a boson star model, in which two scalar fields with the $2\times 2$ coupling matrix $C_{ij}$. The field equations for a spherical configuration derived in the previous section become as follows:

[TABLE]

where we should remember that $C_{12}=C_{21}$. It is assumed that $C_{11}C_{22}-C_{12}^{2}>0$, $C_{11}>0$, for the positive definite $V(|\phi_{i}|^{2})$.

In the large coupling limit, the exponential asymptotic distribution of the scalar density is suppressed colpi . Thus, for the multi-scalar case, it is notable that there are regions where densities of some species of scalar fields vanish. Assuming the normalized densities $\rho_{1}\neq 0$ for the boson 1 and $\rho_{2}\neq 0$ for the boson 2 in the core region of the boson star (in the vicinity of the coordinate origin), the algebraic equations (31) and (32) become

[TABLE]

and quite equivalently

[TABLE]

Thus, we define the total ‘density’ $S_{2}$ and it can be expressed as

[TABLE]

Using the Poisson equation (30), we obtain the differential equation on $S_{2}$ as

[TABLE]

where the positive constant $\omega_{2}$ satisfies

[TABLE]

The solution of the above equation takes the form

[TABLE]

because the center of the boson star at $x=0$ should be nonsingular. Here, the scale factor $A_{2}$ is a positive constant. Then, the gravitational potential is written as

[TABLE]

in this core region.

In the outer region of the star, where $\rho_{1}\neq 0$ and $\rho_{2}=0$, field equations become

[TABLE]

and taking account of the Poisson equation (30), we find

[TABLE]

where the positive constant $\omega_{1}$ satisfies

[TABLE]

Then the solution can be written as

[TABLE]

where $\delta_{1}$ and $A_{1}$ are constants. The gravitational potential then becomes

[TABLE]

in this outer region.

We define the boundary of two regions, where $\rho_{2}\neq 0$ and $\rho_{2}=0$, as $x=x_{2}$ and define the outermost surface of the boson star as $x=x_{1}$, where $\rho_{1}=\rho_{2}=0$. Then, we find that $\rho_{1}\neq 0$ and $\rho_{2}\neq 0$ for $0\leq x\leq x_{2}$, $\rho_{1}\neq 0$ and $\rho_{2}=0$ for $x_{2}\leq x\leq x_{1}$, and $\rho_{1}=\rho_{2}=0$ for $x\geq x_{1}$.

Thus, at $x=x_{2}$, since the total density is continuous,

[TABLE]

is satisfied.444Note that there is no condition on the first derivative of $\rho_{1}(x)$ and $\rho_{2}(x)$ at $x=x_{2}$. Further, since the gravitational force which is derived from the derivative of the potential $\Phi^{\prime}(x)\equiv\frac{d\Phi}{dx}$ varies continuously even at the boundary, the equality

[TABLE]

should be hold. Combining two equalities (45) and (46), we obtain

[TABLE]

and this equation tells us the value for $\delta_{1}$ if $x_{2}$ is given. On the other hand, the outer radius of the boson star $x_{1}$ is determined by the simple equation

[TABLE]

At last, the analytic solution can be obtained, for a given $x_{2}$, as follows.

[TABLE]

[TABLE]

where $\delta_{1}$ is considered as the solution of Eq. (47), also hereafter in the present section. Note that $S(0)=\rho_{1}(0)+\rho_{2}(0)$.

The Newtonian gravitational potential is solved as

[TABLE]

where we set $\Phi(\infty)=0$ by using the shift invariance (15) for potentials.

The chemical potentials can also be obtained as

[TABLE]

which are is always negative, as for a bound state.

A typical profile of the bi-scalar boson star and the gravitational potential are exhibited in FIG. 2 and FIG. 2, respectively, where the couplings are set as $C_{11}=C_{22}=1$, $C_{12}=0.9$ and $x_{2}$ is chosen as $x_{2}=0.35$. In FIG. 2, one can find that the density profile has a kink structure at $x=x_{2}$. This is due to the approximation of the large coupling limit, since the approximation is equivalent to assuming strong but very short-range repulsive forces among bosonic particles. As we will see later, density profiles seem to be almost smooth in multi-scalar cases.

Due to the scale invariance, we obtain a general solution by multiplying an arbitrary common constant with the above set of the solution for $\rho_{1}(x)$, $\rho_{2}(x)$, $u_{1}$, $u_{2}$ and $\Phi(x)$.

The fraction of the particle number of the boson 2, which lives inner region of the boson star, is expressed as

[TABLE]

where $n=n_{1}+n_{2}$. Note that, for a boson star solution, $\omega_{2}x_{2}<\pi$ should hold.

The ratio of the Newtonian binding energy to the total mass of the boson star can be expressed as

[TABLE]

This expression shows the binding energy $E$ is always negative for a boson star solution.

We show the fraction of the particle number of the boson 2 as the function of $x_{2}$ in FIG. 5 and the ratio of the Newtonian binding energy to the total mass of the boson star as the function of $x_{2}$ in FIG. 5, when $C_{11}=C_{22}=1$, $C_{12}=0.9$. In FIG. 5, we show the ratio of the Newtonian binding energy to the total mass of the boson star as a function of the fraction of the particle number of the boson 2, in the same case. We find that the ratio of the Newtonian binding energy to the total mass of the boson star is nearly constant for $n_{2}/n>0.1$. Note that $n_{2}/n$ is independent of the overall scale factor while $E/(mn)$ is proportional to the overall scale.

Now, we consider the gravitational potential and the circular velocity. When we vary the value of $x_{2}$, the shape of the boson star varies and the gravitational potential varies at the same time. Because of the scale invariance under (21), we should focus our mind on the shape, not on the amount.

In FIG. 7, we show the various profiles of the boson stars for $x_{2}=0.3,0.35,0.4,0.45,0.5$, when $C_{11}=C_{22}=1$, $C_{12}=0.9$. The rotation speed $v$ of any object in the circular motion with radius $x$ in the Newton potential $\Phi(x)$ is proportional to $\sqrt{x\Phi^{\prime}(x)}$. In the vacuum region, $\Phi(x)$ is proportional to $-M/x$, where $M$ is the total mass of the boson star. Then, the rotation velocity $v$ becomes $\propto\sqrt{M/x}$ outside the boson star. We exhibit $v/\sqrt{S(0)}=\sqrt{x\Phi^{\prime}(x)}/\sqrt{S(0)}$ for $x_{2}=0.3,0.35,0.4,0.45,0.5$ in FIG. 7. Because the potential spreads out by the density tail due to $\rho_{1}(x)$, the rotation curve can have a flat region, especially around $x_{2}\approx 0.45$ in this case. If a single scalar field model is considered, which is realized by $\rho_{2}(x)\equiv 0$, the range of the gravitational potential becomes narrow and the rotation curve looks far from a satisfactory explanation of the observational data. It is interesting to point out that if the fraction of boson 2 is much larger ($\approx 0.5$, when $x_{2}\approx\pi\omega_{2}^{-1}$), the profile of the boson star becomes much alike a single-scalar boson star. This fact can be read from the rotation curve in FIG. 7.

So far, we have picked up an example of the couplings $C_{11}=C_{22}=1$ and $C_{12}=0.9$. In this case, $\omega_{2}=6.32456$ and $\omega_{1}=1.41421$. The density profile near $\rho_{1}\approx 0$ decreases almost linearly in $\omega_{1}x$. On the other hand, the width of distribution in the core region is determined by $\omega_{2}^{-1}$. Therefore, the broad gravitational potential can be obtained if $\omega_{2}/\omega_{1}\gg 1$. In the above example, we find $\omega_{2}/\omega_{1}=4.472$.

To see this necessary condition more closely, we take another example with the couplings $C_{11}=C_{22}=1$ and $C_{12}=0$. In this case, we find $\omega_{2}=2$ and $\omega_{1}=1.41421$. The density profiles and the rotation curves for various values for $x_{2}$ ($x_{2}=0.8\sim 1.6$) are exhibited in FIG. 9 and FIG. 9, respectively. From this result, we confirm the necessity of the ‘hierarchy’ in the inverse length scales $\omega_{1}$ and $\omega_{2}$.

In this section, we consider the simplest bi-scalar model. We can find that the gravitational potential is spread out by the existence of $\rho_{1}$, and the tail of the boson star leads to an improvement of rotation curve of the galaxy.555As is known, the flat curve is realized when $S(x)\propto x^{-2}$. It is however difficult to obtain the realistic galaxy rotation curve being fit to the observational data by manipulating such a simple model. In the next and subsequent sections, we examine multi-scalar models, since the rotation curve may have the multi-scale structure.666We consider that the density of the scalar field is the total galactic mass density, as the zero-th order approximation, in this paper. The effect of the gravitational mass other than the boson star is briefly discussed in Appendix A.

V Spherical boson stars of multi-scalar fields with a diagonal coupling

matrix

Hereafter, we consider self-interacting $[U(1)]^{N}$ scalar field theory. In this section, we examine the case that $[U(1)]^{N}$ symmetric scalar fields have only self-interactions in individual fields, and no mutual interaction to other scalar fields. Simply speaking, we consider the case with the coupling matrix expressed as a diagonal matrix in the present section, i.e., $C_{ij}=0$ if $i\not=j$.

To obtain the spherical static boson star solution, the equations we should solve are now

[TABLE]

[TABLE]

Thus, $\rho_{i}$ is given by

[TABLE]

We take $\rho_{i}\neq 0$ for $i=1,\dots,k$ and $\rho_{k+1}=\cdots=\rho_{N}=0$ $(1\leq k\leq N)$ without loss of generality. We then define

[TABLE]

Operating the Laplacian on the both sides of the above equation and using the Poisson equation, we obtain

[TABLE]

We assume that $\rho_{i}(x)$ vanishes at $x=x_{i}$ and the outermost surface of the boson star is located at $x_{1}$, where $0=x_{N+1}<x_{N}<x_{N-1}<\dots<x_{2}<x_{1}$. All the normalized density $\rho_{i}$ take nonzero values in the region $0\leq x\leq x_{N}$. In the region $x_{k+1}\leq x\leq x_{k}$, $\rho_{N}=\dots=\rho_{k+1}=0$ and then $S(x)=S_{k}(x)$.

The general solutions of Eq. (60) are

[TABLE]

where $\delta_{k}$ and $A_{k}$ are constants. Note that $\delta_{N}=0$ because of the regularity at the origin.

The gravitational potential is then given by

[TABLE]

Because of the continuity of $S(x)$ and $\Phi^{\prime}(x)$, we find the condition

[TABLE]

This is the recursive equation to determine the value of $\delta_{i-1}$ from $\delta_{i}$ when the other parameters are given.

First of all, we take the simplest case, $C_{ii}=1$ and $C_{ij}=0$ for $i\not=j$, i.e., the coupling matrix $\mathcal{C}=(C_{ij})$ is the identity matrix. In this case, because

[TABLE]

the relation $\omega_{N}>\cdots>\omega_{1}$ holds. As we have seen in the previous section, we need a ‘hierarchy’ in $\omega$ to obtain the boson star profile with a small core and a long tail. For sufficiently large $N$, this condition is satisfied since $\omega_{N}/\omega_{1}=\sqrt{N}\gg 1$.

We show a typical case of $N=10$ in FIG. 11 and FIG. 11, where $\{x_{10},x_{9},x_{8},x_{7},x_{6},x_{5},x_{4},x_{3},x_{2}\}=\{0.43,0.48,0.53,0.58,0.63,0.68,0.73,0.78,0.83\}$. In this case, $\omega_{10}/\omega_{1}=3.16$.

For the sake of comparison, we show the case with $N=5$ in FIG. 13 and FIG. 13, where $\{x_{5},x_{4},x_{3},x_{2}\}=\{0.85,0.9,0.95,1\}$. Because $\omega_{5}/\omega_{1}$ equals $1.73$, which is not so large enough, a small core and a long tail can hardly be obtained even if we choose a smaller value for $x_{5}$, so on.

Next, we consider the case that the coupling matrix is not the identity matrix but a general diagonal matrix. As for a typical case, the profile of a five-scalar boson star and the rotation curves are shown FIG. 15 and FIG. 15, respectively. Here, we take $\mathcal{C}=diag.(5,4,3,2,1)$ and $\{x_{5},x_{4},x_{3},x_{2}\}=\{1.25,1.35,1.45,1.55\}$. Then we find that $\omega_{5}/\omega_{1}=3.38$.

In this section, the large number of scalar fields and/or the large hierarchy in the couplings is necessary for a flat rotation curve, in the case of the diagonal coupling matrix. Both the large number of the fields and the hierarchy in couplings ruin the simpleness of the original theory, which is implicitly desired in theoretical models.

In the next section, we investigate more general coupling matrix including mutual interactions among scalar fields. Since we prefer the model with some symmetry, we propose models associated with a certain symmetric ‘graph’, which appears in graph theory, in the next section.

VI a mutually-coupled scalar model associated with a graph Laplacian

We start with a general case with non-diagonal coupling matrix $\mathcal{C}=(C_{ij})$. We assume that $\rho_{i}>0$ for $i=1,\dots,k$ and $\rho_{k+1}=\dots=\rho_{N}=0$ in the region $x_{k+1}<x<x_{k}$. Then the algebraic equation for nonzero $\rho_{i}(x)$ can be expressed as

[TABLE]

where $\mathcal{C}_{[k]}$ is a $k\times k$ principal submatrix of the matrix $\mathcal{C}$ defined as

[TABLE]

whereas $\mathbf{w}_{[k]}$, $\mathbf{\rho}_{[k]}(x)$ and $\mathbf{u}_{[k]}$ are vectors with $k$ components given by

[TABLE]

Note that $\mathcal{C}_{[N]}=\mathcal{C}$.

In the region $x_{k+1}<x<x_{k}$, $\mathbf{\rho}_{[k]}(x)$ can be solved as

[TABLE]

where $\mathcal{C}_{[k]}^{-1}$ is the inverse of $\mathcal{C}_{[k]}$. Then, we find

[TABLE]

As in the previous sections, by using the gravitational Poisson equation, we obtain

[TABLE]

where

[TABLE]

Now, since we obtain $\omega_{k}^{2}$, we can evaluate the profile of the boson star and the gravitational potential by the same method as in the previous sections.

In the rest of this section, we consider a concrete model whose coupling matrix is represented by a graph Laplacian, which appears in texts of spectral graph theory mohar1 ; mohar2 ; mohar3 ; merris ; GR ; CRS .

Let $G(\mathcal{V},\mathcal{E})$ be a graph with a vertex set $\mathcal{V}$ and an edge set $\mathcal{E}$. The set of edges connects the vertices. A pair of vertices $v$ and $v^{\prime}$ are said to be adjacent, denoted $v\sim v^{\prime}$, if there exists an edge which connects $v$ and $v^{\prime}$. The degree of a vertex $v$, denoted $deg(v)$, is the number of edges directly connected to $v$.

The graph Laplacian $\triangle(G)$ of the graph $G$ is defined by

[TABLE]

where $v,v^{\prime}\in\mathcal{V}$.

Now, we take a model whose coupling matrix can be written as

[TABLE]

where $\mathcal{I}$ is the $N\times N$ identity matrix and $\gamma$ is a non-negative constant. Here, we suppose that the graph $G$ has $N$ vertices. We inversely find the scalar quartic interaction in this model has the form

[TABLE]

where the second sum is done once over all adjacent pairs, and then we find the potential is obviously positive-semidefinite. The use of the graph Laplacian guarantees the positivity of the potential energy in the model.

First of all, we adopt a complete graph $K_{N}$ as a maximally symmetric graph. In a complete graph, each vertex is adjacent to all the other vertices (FIG. 16). The model based on a complete graph is most democratic, because the potential is invariant under any exchange of scalar field species, in other words, vertices in a complete graph have a symmetry under the symmetric group $S_{N}$. The graph Laplacian of the complete graph with $N$ vertices is

[TABLE]

To obtain $\omega_{k}^{2}$, we have to calculate the matrix $\mathcal{C}_{[k]}^{-1}$. We consider the eigenvectors of the matrix $\mathcal{C}_{[k]}$, which is denoted by $\mathbf{v}_{[k]}^{a}$ $(a=0,\dots,k-1)$. They satisfy

[TABLE]

where $\lambda_{[k]}^{a}$ is the eigenvalues of the matrix $\mathcal{C}_{[k]}$. Then, we can express the inverse matrix using the eigenvalues and the eigenvectors of the matrix as

[TABLE]

In the present case that $\mathcal{C}=\mathcal{I}+\gamma\triangle(K_{N})$, the vector $\frac{1}{\sqrt{k}}\mathbf{w}_{[k]}$ is an eigenvector of $\mathcal{C}_{[k]}$ for any $k$, and the associated eigenvalue is $1+\gamma(N-k)$. Therefore,

[TABLE]

where we used the orthogonal relation among the eigenvectors.

The necessary condition for a ‘small core and long tail’ boson star solution is that the ratio $\omega_{N}/\omega_{1}$ is much larger than unity, actually $\omega_{N}/\omega_{1}\stackrel{{\scriptstyle>}}{{\sim}}3\sim 4$. In the present model, we find

[TABLE]

so, the choice of $\gamma\simeq 1$ and $N\simeq 3\sim 4$ satisfies the condition.

We demonstrate the calculation for $N=5$, $\gamma=1$, and $\{x_{5},x_{4},x_{3},x_{2}\}=\{0.83,1.03,1.23,1.43\}$, and show the results in FIG. 18 and FIG. 18. By integrating densities, we can evaluate the particle number of each scalar boson. The fraction of species are $\{n_{1}/n,n_{2}/n,n_{3}/n,n_{4}/n,n_{5}/n\}=\{0.786,0.0788,0.0620,0.0453,0.0283\}$, where $n=\sum_{i=1}^{5}n_{i}$. Although the value of $n_{1}/n$ is larger than those of others, because of the long tail distribution of the boson star, this composition is not so unnatural. The contribution of each boson species are shown in FIG. 19.

Thus, we have found that the model with the coupling matrix associated with the graph Laplacian of the complete graph leads to flat rotation curves in a natural setting, for comparatively small number of scalar fields which have symmetry under any exchange of species. Moreover, a small but finite number of fields brings about a diversity of galaxy rotation curves, as observations indicate. It can be recognized that the variation on diverse galaxies is due to the set of fractions of the particle numbers.

Next, we consider the use of the other graphs in the coupling matrix. We now consider a cycle graph $C_{N}$ (FIG. 21). The coupling matrix is given in the same form Eq. (73), as previously considered. The model based on a cycle graph still has discrete symmetry under $\phi_{i}\rightarrow\phi_{i+l}~{}(\mbox{mod }N)$, where $l$ is an integer.

The graph Laplacian of a cycle graph $C_{N}$ takes the form

[TABLE]

In the model associated with $\triangle(C_{N})$, we find

[TABLE]

Note that $\frac{1}{\sqrt{N}}\mathbf{w}_{[N]}$ is the eigenvector belonging to zero eigenvalue for any simple graph Laplacian. Therefore, we find777If $deg(v)=d$ for all $v\in\mathcal{V}$, the graph is called a $d$-regular graph. The cycle graph is a $2$-regular graph. For the model associated with $d$-regular graph, one can find $\omega_{N}/\omega_{1}=\sqrt{N(1+\gamma d)}$.

[TABLE]

In general, the necessary condition in this model is more severe than that in the model associated with a complete graph, but a sufficiently large number of scalar fields satisfies the condition. In a loosely connected graph such as a cycle graph, each degree of the vertex is much smaller than $N-1$. Thus, $\omega_{1}$ tends to be larger than the case with a complete graph. In the star graph with $N$ vertices (FIG. 21), the maximal degree of the vertex is $N-1$, which is the same as that in a complete graph. Then, in the model based on the star graph, the ratio ${\omega_{N}}/{\omega_{1}}$ takes the same value as that in the model of the complete graph $K_{N}$ with the same number $N$ of scalar species. However, the model associated with the star graph with $N$ vertices has less symmetry, i.e. symmetric under the symmetric group $S_{N-1}$ rather than $S_{N}$.

Recently, models with mass hirarchy generated by the ‘clockwork mechanism’ CI ; KaRa ; GM1 ; FPRT ; KeRi ; HTT ; teresi ; AD ; CGS ; GM2 are eagerly investigated by many authors. If it is feasible to use the same structure in our coupling matrix instead of their mass matrices, we would have an interesting model for boson stars constructed by several fields.

VII Summary and outlook

In the present paper, we have examined the Newtonian boson star with many $U(1)$ charges in the large coupling limit. The explanation of rotation curves of galaxies by gigantic boson stars is improved in this model.

The necessary condition for a flat rotation curve is that $\omega_{N}$ is several times larger than $\omega_{1}$ in terms of our model parameter. This is naturally led from the model with a large number of scalar field and/or the model whose coupling matrix is associated with the graph Laplacian of the complete graph. The latter model has larger symmetry on several scalar fields. It is worth pointing out that several species of bosons are necessary to make a variety in the rotation curves of different galaxies, even including dwarf galaxies.888We may need, however, a large number $N$ naturally to obtain hierarchical scales.

In the present paper, we have concentrated ourselves on the galactic boson stars. The structure of small multi-scalar boson stars is still similar because of the scale invariance of our models. It will be interesting to find some new aspects of the boson stars in various circumstances, such as in a collision of multi-scalar boson stars.999The collision of boson stars has been studied in Refs. BG ; GG , for example.

Although the model we have demonstrated is rather a toy model and has only been focused as in the large coupling limit, it is the simplest effective theory of multi-scalar boson stars. As a variant of the model, it is interesting to consider the case that there exist other fields which do not participate the ingredient of the boson star but affect on the creation or decay of boson stars. Anyway, we would like to find the relevance of the multi-scalar theory to particle physics, and wish to clarify its role or significance.

In future work, we will consider general relativistic boson stars and graph-oriented models with many charges or arbitrary couplings. We also have much interest on boson stars in plausible models in which finite or zero couplings. It is known that single-scalar boson stars with finite self-couplings can also be approximated analytically by connecting the exponential tail in the asymptotic region of the boson star GA . Since the rotation curve is most sensitive on the tail of the boson star, analyses of models with finite self-couplings should be performed as in the next step. Time dependent solutions or oscillations of multiple scalar fields are also interesting. We wish to investigate these subjects elsewhere.

Appendix A a multi-scalar boson star with an external gravitational source

The realistic galaxies have bulges and halos of stars and gas, which are gravitational sources and are expected not to interact with scalar fields in our model (and other unknown dark stuff).

In the presence of the nonnegligible external gravitational source,101010We neglect further back reaction to the source from the gravitation of scalar bosons. whose energy density is given by $\rho_{ext}$ other than the scalar fields, we add an inhomogenious term in the right hand side of Eq. (56) as

[TABLE]

where $\bar{\rho}_{ext}=4\pi G\Lambda/m^{2}\rho_{ext}$. Although the normalization seems to be strange, the ratio of energy densities are nevertheless simply given by $\bar{\rho}_{ext}/S$. Now, we should solve the differential equations

[TABLE]

and the special solution for the inhomogeneous term is found to be

[TABLE]

Since a main subject of the present paper is an interest in analytical study, we consider a simple example in this appendix. We now assume

[TABLE]

where $A_{e}$ and $\omega_{e}$ are constants and $\bar{\rho}_{ext}(x)=0$ for $x>\pi/\omega_{e}\equiv x_{N+1}$. We further assume $\omega_{e}>\omega_{N}$ for simplicity. Then, we find

[TABLE]

where $A_{N}$ is an arbitrary constant. We have only to find the total solutions by the connection condition of $-x\Phi^{\prime}/S$ as in Sec. V.

A result in the model associated with $K_{5}$ in Sec. VI is shown in FIG. 23 and FIG. 23, where we set $\omega_{e}=\sqrt{5}\omega_{5}$, $\bar{\rho}_{ext}(0)/S(0)=0.5$ and $\{x_{5},x_{4},x_{3},x_{2}\}=\{0.44,0.83,0.95,1.07,1.19,1.31\}$. The most appropriate parameter set for a flat rotation curve, of course, changes slightly by the effect of the gravitational source. Conversely, it would be said that the density profle from the present model can be adjusted to various distributions of ordinary matter in many galaxies.