{'fraction_non_alphanumeric': 0.03666092423006014, 'fraction_numerical': 0.012149786925378547, 'mean_word_length': 4.392438070404172, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Numerical studies of Phi^2-Oscillatons

We present an exhaustive analysis of the numerical evolution of the Einstein-Klein-Gordon equations for the case of a real scalar field endowed with a quadratic self-interaction potential. The self-gravitating equilibrium configurations are called oscillatons and are close relatives of boson stars, their complex counterparts. Unlike boson stars, for which the oscillations of the two components of the complex scalar field are such that the spacetime geometry remains static, oscillatons give rise to a geometry that is time-dependent and oscillatory in nature. However, they can still be classified into stable (S-branch) and unstable (U-branch) cases. We have found that S-oscillatons are indeed stable configurations under small perturbations and typically migrate to other S-profiles when perturbed strongly. On the other hand, U-oscillatons are intrinsically unstable: they migrate to the S-branch if their mass is decreased and collapse to black holes if their mass is increased even by a small amount. The S-oscillatons can also be made to collapse to black holes if enough mass is added to them, but such collapse can be efficiently prevented by the gravitational cooling mechanism in the case of diluted oscillatons.


I. INTRODUCTION
In a seminal paper Seidel & Suen [1] found that there exist non-singular, time-dependent equilibrium configurations of self-gravitating real scalar fields. These oscillating soliton stars are called oscillatons, and are solutions of the Einstein-Klein-Gordon (EKG) system of equations for minimally coupled real scalar fields. The time-dependence of these solutions appears as a fundamental ingredient that allows singularities to be avoided, in contrast to static solutions of the EKG equations with real scalar fields, where singularities frequently appear. The case of real scalar fields is quite different to the case of complex scalar fields, for which the EKG equations provide the so-called boson stars, which are non-singular solitonic solution with a static geometry (the components of the complex field oscillate in precisely such a way that the stress-energy tensor is time independent). However, the stability of oscillatons has not been studied in as much detail as in the case of boson stars [2,3,4]. Such studies are necessary because of the possible role of oscillatons in astrophysics and cosmology, where real scalar fields have been proposed as candidates for the dark matter in the Universe [5,6,7].
of the evolution of the EKG equations, much in the same way as it has been done for boson stars [2,3,4]. The main aim is to investigate whether oscillatons are stable. For simplicity, we only consider here the spherically symmetric case of a real scalar field Φ endowed with a quadratic scalar potential of the form V (Φ) = (m 2 /2) Φ 2 . Other cases will be treated in future publications. This paper is organized as follows. In section II, we present the necessary mathematical background to find equilibrium configurations and to further evolve the EKG equations. As we shall see, the equilibrium configurations may be classified in two different groups: the S and U-branches. The numerical methods and tests are presented in section III. In sections IV and V, we analyze the results obtained from the evolution of the EKG equations for different oscillatons. We separately study the behavior of S and U-oscillatons by adding to them small and strong perturbations that change their total mass. We conclude in section VI.

II. MATHEMATICAL BACKGROUND
To begin with, we consider the spherically symmetric line element ds 2 = −α 2 dt 2 + a 2 dr 2 + r 2 dθ 2 + sin 2 (θ)dϕ 2 (1) with α(r, t) the lapse function and a(r, t) the radial metric function. We choose the polar-areal slicing condition (i.e. we force the line element to have the above form at all times, so that the area of a sphere with r = R is always equal to 4πR 2 ); this choice of gauge will force the lapse function α(r, t) to satisfy an ordinary differential equation in r. Through out the work, we will be using units such that c = = 1, and we express the gravitational constant in terms of he Planck mass: G = 1/m 2 Pl . The energy momentum tensor of a scalar field Φ endowed with a quadratic self-interaction potential V (Φ) = (m 2 /2)Φ 2 is The non-vanishing components of T µν are and also T 3 3 = T 2 2 . These different components are identified as the energy density ρ Φ , the momentum density P Φ , the radial pressure p r , and the angular pressure p ⊥ . The parameter m is interpreted as the mass of the scalar particles. Over-dots denote ∂/∂t and primes denote ∂/∂r. The evolution of the metric functions α and a can be obtained from the Einstein equations G µν = κ 0 T µν , with κ 0 = 8πG. In order to write appropriate evolution equations, we now introduce the first order variables Ψ = Φ ,r and Π = aΦ ,t /α. Also, we define the dimensionless quantities r = x/m, t → t/m, Φ → Φ/ √ κ 0 , where we note that the bosonic mass is the natural scale for time and distance.
Using these new variables, the Hamiltonian constraint becomes and the polar-areal slicing condition takes the form The evolution of the scalar field Φ is governed by the Klein-Gordon (KG) equation, which appears as a consequence of the conservation equations of the scalar field energy-momentum tensor (2) in the form where ✷ = g αβ ∇ α ∇ β is the d'Alambertian operator. The KG equation is equivalent to the following set of first order differential equations Equations (7)(8)(9)(10)(11)(12) form the complete set of differential equations to be solved numerically. The evolution equation for Π above is further transformed into the equivalent form: Notice that the first term on the right hand side of this equation includes now a first derivative with respect to x 3 (and not a third derivative). The reason for doing this transformation has to do with the numerical regularization near the origin of the 1/x 2 factor in equation (11) above (see Ref. [11]).
In order to deal with the non-linearities present in the EKG equations, it is convenient to introduce the new variables A(t, r) = a 2 (t, r), C(t, r) = (a/α) 2 , for which eqs. (7)(8)(9)(10)(11)(12) take the form The lapse function is later obtained as α 2 (t, x) = A(t, x)/C(t, x). Notice that Eq. (17) is a consequence of the momentum constraint (the {t, r} part of the Einstein equations.) We shall consider the Fourier expansions where ω is called the fundamental frequency and j max is the mode at which the Fourier series are truncated. Solutions are obtained by introducing the Fourier expansions (18) in eqs. (14)(15)(16)(17), and setting each Fourier coefficient to zero; that is, the EKG equations are reduced to a set of coupled ordinary differential equations.
The boundary conditions are determined by requiring non-singular and asymptotically flat solutions, for which equations (14-17) become an eigenvalue problem. Thus, it is only necessary to determine the initial values φ i (0), C i (0) (the fundamental frequency is an output value) corresponding to a given central central value φ 1 (0), to obtain different equilibrium configurations.  18)). The convergence of the Fourier series is manifest (notice the re-scaling of the higher Fourier modes).
A typical oscillaton solution is shown in Fig. 1. Even though we are solving non-linear equations, the Fourier series converges rapidly. A particular feature of the solutions is that they are represented only by odd Fourier coefficients of the scalar field Φ and the even coefficients of the metric functions A and C.
In Fig. 2, we show the calculated total mass (M T ), the fundamental frequency (Ω ≡ ω/m) and the radius at which the radial metric coefficient reaches its maximum value at t = 0 (R max (0)) for different configurations. In the case of oscillatons, the position of the maximum of a 2 is not a fixed value but instead oscillates in time. However, as we shall see below (see Fig. 8), the amplitude of such oscillations is quite small and the initial value can be taken as representative of each oscillaton. The maximum mass M c = 0.607 m 2 Pl /m Φ is reached for a central value φ 1 c (0) = 0.48, to which also corresponds a fundamental frequency Ω = 0.864. The fundamental frequency is always such that Ω ≤ 1 for all oscillatons, with Ω = 1 for the trivial solution. In general, more massive oscillatons oscillate with a smaller fundamental frequency.

III. NUMERICS
A. Numerical algorithm for the evolution of the system In order to integrate the Klein-Gordon equations (10,11,12) we used a method of lines with second order centered differences in space. For the time integration we use a method inspired in the Iterative Crank Nicholson (ICN) scheme with three iterations (see for example [12,13]). The standard ICN is used to integrate a system of evolution equations of the form: Given values of the variables u i at time step t = n∆t, one updates the values of the variables to time step t = (n + 1)∆t in the following way: where u n i = u i (t = n∆t), u n+1 i = u i (t = (n + 1)∆t), and u (k) i are intermediate values. Taking N = 3 is enough to obtain a second order accurate, stable scheme [12,13] (in fact, taking N = 2 is enough for second order accuracy, but it is unstable).
For our purposes we have modified the above algorithm for the case N = 3 in the following way: The reason for this modification is that the above scheme is considerably less dissipative than standard three-step ICN. In fact, for linear equations one can show that the above scheme is third order in time. We have found our modified scheme ICN to be very stable and robust in practice.
Once we have advanced the variables Φ, Ψ and Π one time step using the above algorithm, we substitute their new values into equations (7) and (8). These are simple ordinary differential equations on the radial coordinate that we solve using a standard second order Runge-Kutta scheme.

B. Boundary conditions
Our set of equations is singular at x = 0. To avoid the singularity, we stagger the origin and take a spatial grid of the form x i = (i − 1/2)∆x. The fictitious point at x 0 = −∆x/2 is used to impose appropriate parity conditions: Π is even and Ψ is odd. Notice that we can integrate Φ all the way to the boundary point x 0 = −∆x/2, since its evolution equation does not require the evaluation of spatial derivatives.
At the outer boundary we also need to impose boundary conditions. Notice again that we do not need to apply a boundary condition for Φ, as its evolution equation can be integrated all the way to the boundary point. For Π we assume that, for large enough x, it behaves as an outgoing wave pulse of the form: with u an arbitrary function. In differential form this becomes which, performing finite difference can be solved for the unknown boundary value at the new time level. On the other hand, it is not difficult to convince oneself that the function Ψ does not behave as an outgoing wave at the boundary. However, we are assuming that Π does and, as a consequence, so does Φ. In particular, the outgoing wave boundary condition applied to Φ can be seen to imply that at the boundary: This equation can then be used to obtain boundary values for Ψ once those of Φ and Π are known. Finally, we need to mention the boundary conditions used for the ordinary differential equations that have to be solved to find the metric functions α and a. For a, we use the fact that local flatness implies that a(x = 0) = 1 and ∂ x a(x = 0) = 0. These two conditions imply that a(x 0 ) = a(x 1 ) = 1+O(∆x) 3 . We use these two boundary values at the first two grid points and integrate the second order hamiltonian constraint outwards.
For the lapse function α we use the fact that, on a vacuum, our slicing condition implies that we are in Schwarzschild coordinates, so we must have α = 1/a. We then assume that our boundaries are sufficiently far away as to be always in a vacuum, and impose α = 1/a as an outer boundary condition. The slicing condition is then integrated inwards. One could presumably improve on this by setting α(x = 0) = constant, ∂ x α(x = 0) = 0 as boundary conditions on the origin, integrating the slicing condition outwards, and then re-scaling it so that the lapse goes as 1 + k/x far away (k is a constant; notice that the slicing condition is scale invariant).

C. Code tests
In order to illustrate how all these ingredients work together properly, we study now the accuracy of our numerical methods in some particular cases. In Figure 3 we show the convergence of a system initially consisting of a Gaussian pulse of scalar field. The equation we use to calibrate our techniques is the (r, t) component of the Einstein's equations (the momentum constraint): which should be satisfied for an exact solution [14]. What we show in the plot is the L2 norm of the value of β across the grid as a function of time for three different resolutions. The initial Gaussian corresponds to an oscillaton with M = 0.575, so that it is in the stable branch and has a long time life; the boundary is located at x = 50 so that the evolutions were carried out up to 100 crossing times. The fact that the value of β goes down by a factor of 4 every time the resolution is doubled shows that the code remains second order convergent.

IV. EVOLUTION OF OSCILLATONS: THE S-BRANCH AND THE QUASI-NORMAL MODES
We shall call S-branch oscillatons those equilibrium configurations to the left hand side of the critical configuration in the plot M vs φ 1 (0) in Fig. 2 (or equivalently, those configurations to the right hand side in the plot of M vs R max ). As a first step in the study of the stability of oscillatons, we will numerically evolve these equilibrium configurations.
We should keep in mind that, compared to realistic oscillatons, the solutions presented in Fig. 2 are already perturbed. These perturbations arise because of the inherent discretization error of the numerical solutions and (mainly) of the truncated Fourier expansions. Therefore, the eigen-solutions of section II A will be treated as already slightly perturbed profiles [1].
The main results presented in this section can be summarized as follows. i) S-oscillatons are stable against small perturbations. ii) Besides the fundamental oscillations of the system, there is an overall vibration of the oscillaton when slightly perturbed. The frequencies (f ) of these vibrations should be identified with the so-called quasi-normal modes of the system, and the resulting plot f vs M is a particular feature of the Φ 2 -oscillatons. As we shall see below, such quasi-normal modes are very important in the study of evolved profiles.
We take the scalar field profiles obtained by solving eqs. (17-16) to be the initial data for the evolution equations. In this form, the initial conditions for the scalar field are (see eqs. (18)) and then the metric functions α, a are calculated through the Einstein equations. In this manner, we will also check the consistency between the eigenvalue problem and the numerical evolution since no additional information (like the metric functions, the fundamental frequency ω, etc.) is taken initially. As a typical example of a S-branch oscillaton, we show in Fig. 4 the evolution of the initial profile shown in Fig. 1, characterized by the central value φ 1 (0) = 0.2828. From the figure we see that such an oscillaton is stable, maintaining the same oscillatory pattern from t ≃ 300 (see Fig. 5) up to times t ≃ 20, 000 (in the plot we are showing only a fraction of the run). That the evolution is stable can also be seen from Fig. 5, in which we show the accuracy of the numerical evolution by plotting the momentum constraint. In Fig 6 we show the evolution of the total integrated mass. We can see a small appreciable adjustment of the original mass at around t = 300, indicating a small ejection of scalar field from the system. The overall linear decay of the mass can be shown to be consistent (using convergence tests) with a small amount of numerical dissipation still present in our numerical method, and is therefore not an intrinsic decay of the oscillaton.
The simulations show that the evolved profile oscillates periodically at two different time scales. In particular, there is a short-period oscillation, which we may identify with the fundamental oscillation, and an overall vibration with a longer period [1]. Following the literature of boson stars [2], we shall refer to the later as the "quasi-normal modes".
To determine the oscillatory scales of the solution, we calculated the power spectrum of the evolution of (g rr ) max in Fig. 4, and plotted it in Fig. 7. Notice that the dominant frequencies are the quasi-normal one and  ). The initial configuration is that of Fig. 1 with ∆x = 0.01, ∆t/∆x = 0.5. The boundary is at x = 100 and the evolution is shown to t = 5, 000 m −1 , some 50 crossing times. In fact, we have followed the run up to t = 20, 000 m −1 and no change was observed. Thus, we can conclude that the initial configuration is stable against small perturbations.
(Bottom) The density profile ρΦ(t, x) times x 2 is plotted for late times.  Higher modes corresponding to even-multiples of the fundamental frequency are also present, but their power is smaller. Observe that the fundamental frequency is a little bit readjusted too, but its value is consistent with the eigen-solution in Fig. 1.
We systematically evolve all other S-branch equilibrium configurations shown in Figs. 2, calculating their corresponding power spectrum and then their quasinormal and fundamental frequencies. In all cases, the evolved profiles behaved accordingly with the description given above for the (representative) case φ 1 (0) = 0.2828.
For instance, the values of the fundamental frequencies were just a little bit readjusted in all cases, but its values were always consistent with those shown in Figs. 2. In general, the power spectrum of more dilute oscillatons is dominated by the fundamental frequency, while the quasi-normal frequency dominates for oscillatons near the critical point. Also, the linear numerical dissipation of mass was as in Fig. 5, suggesting also that it should be attributed to the numerical method and not to the intrinsic properties of oscillatons.
The resulting plot f vs M in Fig. 7 has been very useful to analyze the evolution of boson stars [2,4], and we will see it is useful for oscillatons as well. The plot shows a maximum and a sharp decline near the critical mass. This is a typical behavior indicating the transition from stable to unstable configurations.
'To further demonstrate that the S-oscillatons are stable and that the frequencies shown in Fig. 7 are their intrinsic quasi-normal modes, we show in Fig. 8 the evolution of the total mass M T and R max compared to the equilibrium configurations as shown in Fig. 2. It is clear that the slightly perturbed S-oscillatons oscillate with very small amplitudes around the original equilibrium The frequencies of the quasi-normal modes obtained from the evolution of slightly perturbed S-branch oscillatons. We also show the migration of perturbed S-branch and U-branch oscillatons, labelled A-F; the details of their evolution can be found in sections IV 1,V. In the plot, ∆t represents the time intervals at which the evolved mass was measured for each case.
configurations. In other words, they are not migrating to other oscillatons nor decaying. This can be compared to the migrating oscillatons shown in Figs. 10,12,13 below.
The results presented above give evidence for the points (i-ii) outlined at the beginning of this section. We would like to stress here that the important result is that there do exist stable oscillatons (at least in the S-branch). As a side-effect, we have also proven the consistency between the eigenvalue problem and the numerical evolution code.

Perturbed S-oscillatons
Our interest now is to determine whether S-branch oscillatons are stable against strong perturbations. That is, we want to know the conditions under which such oscillatons will either collapse into a black hole, disperse away or form another oscillaton. The latter is a rather interesting possibility since, as we shall show below, it implies the migration of oscillatons.
We have included Gaussian-like perturbations in the original equilibrium configurations, like that shown in Fig. 9, which can be seen just as perturbations to the original mass of the oscillatons. (We shall call original the profiles shown in Figs. 2 and their corresponding parameters). The main results are summarized as follows. iii) If an oscillaton is perturbed in such a way that its mass is less than the critical mass M c ≃ 0.606 (m 2 Pl /m), [15] it migrates to another solution in the S-branch. iv) On the other hand, if the initial mass is larger than M c , the oscillaton can either migrate to an S-oscillaton or collapse to a black hole. The collapse to a black hole can be prevented by the oscillaton ejecting the mass excess, the so-called gravitational cooling mechanism [3]. This mechanism is highly effective for dilute oscillatons. In both cases (iii,iv) above, the evolution of the perturbed profiles can be tracked according to their vibration frequencies, i.e., the quasi-normal modes. We have plotted some typical evolutions of perturbed S-branch oscillatons in Fig. 10 which show the different behaviors described above.
We start by discussing point (iii). In the first case, we decreased the initial mass of a φ 1 (0) = 0.45-oscillaton to an initial value M i = 0.484 (20% less than its original mass M = 0.605). We can observe that initially the oscillaton expands and then bounces back while it looses mass, to finally settle down on the S-branch. The migration is manifest in the trajectory shown in Fig. 10, that oscillates around the stable equilibrium configurations discussed in the previous section. During the evolution, the oscillaton maintains a fixed vibration frequency. Its path of migration is labelled A in Fig. 7, and shows that the oscillaton will stop at a diluted one.
A second case corresponds to a φ 1 (0) = 0.1-oscillaton which was strongly perturbed as to have an initial mass of M i = 0.579 (40% over its original mass M = 0.414). Despite the strong perturbation, the oscillaton collapses and looses enough mass to settle down onto another Soscillaton. By measuring the vibration frequency, we find that the migrating oscillaton follows the path labelled B in Fig. 7, which suggests that the oscillaton is migrating to a φ 1 (0) = 0.15-like profile. This can also be seen in Fig. 11, in which the profile of the metric coefficient g rr rapidly approaches and oscillates around the final configuration.
To illustrate point (iv), a φ 1 (0) = 0.1-oscillaton was perturbed to have an initial mass of M i = 0.687, a value larger than the critical one. However, we can observe that the oscillaton rapidly loses mass and migrates back to the S-branch. Its migration path is labeled C in Fig. 7. This is a typical example of the efficiency of the gravitational cooling mechanism in enough dilute oscillatons.
The opposite case is shown by a perturbed φ 1 (0) = 0.2-oscillaton, whose initial mass is M i = 0.657. The gravitational cooling mechanism is not efficient enough in this case to prevent the formation of a black hole. We see this by noticing that its trajectory stops at the value of R max corresponding to the Schwarzschild radius This mechanism is so efficient in the case of dilute oscillatons, that it can prevent the formation of a black hole even if the initial mass is larger than the critical value. However, this cooling mechanism is sometimes not enough for dense oscillatons. The asterisks denote the starting points. The runs were followed up to t = 5, 000m −1 .
of the hole. However, we have observed that the same oscillaton does migrate to the S-branch if the initial mass is somewhat less perturbed, while still being larger than the critical one. This is also a generic phenomenon we have observed by perturbing other S-oscillatons. This tells us that the gravitational cooling mechanism is highly efficient if the scalar field is diluted enough, like in the cases shown in the original paper of Seidel & Suen [3].
The results presented so far point out that, apart from being stable configurations, S-branch oscillatons are indeed the final states in the evolution of other perturbed S-oscillatons. These two properties (stability and final state-quality) are the imprint of S-oscillatons. To end this section, we would like to stress here the important role of the quasi-normal modes to follow the migration process of oscillatons. The examples presented so far show that the vibration frequency of the system remains the same during its evolution, even while the system looses mass.

V. EVOLUTION OF OSCILLATONS: THE U-BRANCH
We shall call U-branch oscillatons those equilibrium configurations on the right (left) hand side of the critical configuration in the plot M vs φ 1 (0) (M vs R max ) in Fig. 2. We also evolved these equilibrium configurations, under the same idea that the eigen-solutions of section II A are already slightly perturbed configurations. The main result of this section is that U-oscillatons are intrinsically unstable since, as we shall see below, they decay and migrate to the S-branch under small perturbations.
Shown in Fig. 12 are some instances of slightly perturbed U-oscillatons, and the plot M vs R max speaks by itself when compared to Fig. 8. The U-branch oscillatons are intrinsically unstable; they migrate to and settle down onto the S-branch, even under small perturbations. The larger the central value (φ 1 (0) > 0.48), the more unstable they are. Their migration to the S-branch also confirms the stability of the S-oscillatons and their important role as the final states in the evolution of migrating oscillatons.
As we have done before, the migration of oscillatons can be tracked by determining their quasi-normal frequencies. For example, the slightly perturbed φ 1 (0) = 0.7, 0.8-oscillatons shown in Fig. 12 are labeled D,E, respectively, in Fig. 7.

A. Perturbed U-oscillatons
We have also perturbed the U-branch equilibrium configurations and studied whether they migrate or collapse into black holes. The main results of this section are as follows. i) If the initial mass somewhat larger than the original mass of the eigen-configuration, the oscillaton collapses into a black hole. This result is independent of the kind of perturbation made and then it can be characterized by the excess in mass only. ii) If the mass of the original configuration is decreased by the perturbation, the oscillaton migrates to the S-branch.
The φ 1 (0) = 0.8-oscillaton is a typical example of Uoscillaton and its evolution under strong perturbations is plotted in Fig. 13. If the initial mass of this oscilla- ton is perturbed so that is 2% more massive that the original mass, the final stage of the collapse is a black hole, even though the initial mass is less than the critical value. This can also be seen in Fig. 14, where we see plots of the metric coefficients g tt and g rr . The coefficient g tt shows the well-know "collapse of the lapse" (the lapse drops to zero) characteristic of black holes, while the coefficient g rr shows the "grid stretching" effect (te radial metric grows exponentially), also typical of black hole evolutions. On the other hand, if the initial mass of the same oscillaton is 2% less than the original mass, the oscillaton rapidly migrates to the S-branch and settles down into a stable oscillaton. Its evolution path in Fig. 7 (label F), overlaps with that of the migrating S-oscillaton shown in Fig. 9,11 (label B in Fig. 7). Thus, the final state of the evolution seems to be a φ 1 (0) = 0.15-oscillaton, as can also be seen in Fig. 15. This last fact reveals again the usefulness of the plot f vs M to determine the final configuration of a migrating oscillaton. We would like to summarize here the results presented so far. We have found that the imprint of U-branch oscillatons is that they are intrinsically unstable. That is, even small perturbations provoke their migration to other configurations, instead of the small oscillations around the original ones. This latter fact again confirms the stability and final state-status of S-oscillatons. Another particular characteristic of U-oscillatons is that the formation of a black hole is very likely if the original profile is perturbed by adding mass to it. This phenomenon occurs even if the perturbed initial mass is less than the critical one, so that the gravitational cooling mechanism is not efficient at all in the U-branch.

B. The S-U transition point
Now we take a closer look at the equilibrium configurations near the critical profile, the one being the most massive in Figs. 2, which we shall call the S-U transition point.
To begin with, let us recall the basic characteristics of the S,U-branches. The basic imprint of S-branch oscillatons is that they are stable and then oscillate around the initial profiles when slightly perturbed, see Fig. 8. On the other hand, U-branch oscillatons are intrinsically unstable and they migrate to the S-branch even when slightly perturbed, like those cases shown in Fig. 12.
The evolution of the profiles near the S-U transition point when slightly perturbed are shown in Fig. 16. We observe that the results summarized above are again con- firmed, and then the critical configuration characterizes the transition between the S and U branches under small perturbations. We have also found the same results in the case of strong perturbations. If the mass of these oscillatons is decreased, they migrate to S-oscillatons. If the mass is increased even by a small amount, they collapse rapidly forming black holes. That U-oscillatons form black holes is not surprising. The fact that S-oscillatons can also collapse to black holes near the critical point is a consequence of the inefficiency of the gravitational cooling mechanism in the case of dense S-oscillatons, as it was found in section IV 1. Actually, if the mass of an S-oscillaton near the transition point is perturbed to be larger than M > 0.606, it will collapse into black holes; otherwise, it will migrate to more diluted configurations. That is, the value M c ≃ 0.606 can be seen as the true critical value in the S-U transition region.
Thus, we can conclude that the critical oscillaton (φ 1 (0) ≃ 0.48, M c ≃ 0.606m 2 Pl /m) marks the transition from the S-branch to the U-branch.

VI. CONCLUSIONS
In this paper, we have studied the properties of oscillatons by numerically evolving the Einstein-Klein-Gordon equations. According to the results, oscillatons can be classified into two well definite groups, the S and Ubranches.
S-branch oscillatons are stable equilibrium configurations, which oscillate according to their fundamental and intrinsic frequency (ω). When slightly perturbed, these oscillatons vibrate with definite frequencies which we have identified as the frequencies (f ) of their particular quasi-normal modes. This is supported by the fact that these vibrations are small-amplitude oscillations around the equilibrium configurations. If these oscillatons are strongly perturbed such that their original mass is decreased, they migrate to and settle down onto other Sbranch oscillatons. This fact points out that S-oscillatons should be seen as final states to which other scalar configurations migrate to.
On the other hand, we have found that S-oscillatons do not, in general, collapse into black holes if their mass is increased by the perturbation. Even though we have found that oscillatons with high-central densities form black holes if their mass is larger than the critical value M c ≃ 0.606, diluted oscillatons can avoid such fate by means of the gravitational cooling mechanism. The later can be so efficient that strongly perturbed oscillatons can migrate to stable configurations, even if their initial mass was much larger than the critical value.
The evolutions of the U-branch oscillatons consistently show that they are intrinsically unstable configurations. Even under small perturbations, they migrate to the S-branch. This migration also appears when they are strongly perturbed in such a way that their original mass is decreased. However, another manifestation of their unstable nature is that they rapidly collapse into black holes if their mass is increased, even by small amounts and even if the initial mass is well below the critical value.
In all cases of migrating oscillatons, we found that their evolution can be followed by measuring their vibration frequency, and this can be used to predict to which configuration the perturbed oscillatons will migrate to.
Through this paper, we have studied oscillatons with a quadratic scalar potential for simplicity. However, we have already initiated similar studies to include a selfinteraction in the potential and they are work in progress to appear elsewhere.