Oscillation modes of ultralight BEC dark matter cores

TL;DR

This paper investigates the oscillation modes of ultralight bosonic dark matter cores, providing a systematic numerical analysis of their frequencies to better understand their stability and observational signatures in galaxy formation.

Contribution

It introduces a comprehensive numerical study of the oscillation frequencies of equilibrium configurations of ultralight BEC dark matter cores under various perturbations.

Findings

Identified dominant oscillation frequencies for modes l=0 to 4.

Demonstrated stability of configurations under perturbations.

Provided data to link oscillation modes with observable galaxy features.

Abstract

Structure formation simulations of ultralight bosonic dark matter, ruled by the Gross-Pitaevskii-Poisson (GPP) system of equations, show that dark matter clumps to form structures with density profiles consisting of a core surrounded by a power-law distribution. The core has a density profile similar to that of spherically symmetric equilibrium configurations of the GPP system. These configurations have been shown to be stable under a variety of perturbations and to have attractor properties. It is interesting to know the dominant frequencies of oscillation of these configurations. One reason is that in galaxies the effects of perturbations can trigger observable effects for specific frequency modes. Another reason is that during the process of structure formation, the oscillation modes can help to characterize a core. Based on these motivations, in this manuscript we present a systematic numerical analysis of the reaction of equilibrium configurations to various axi-symmetric perturbations with modes $l=0,1,2,3,4$. We then calculate the first few oscillation frequencies of equilibrium configurations for each mode.

Tables

Table 1. Table 1: Frequencies measured along the radial direction with detector D 1 subscript 𝐷 1 D_{1} and along the z − limit-from 𝑧 z- axis with detector D 2 subscript 𝐷 2 D_{2} . Δ Δ \Delta is the standard deviation of the frequency for each peak, calculated with the results of all the combinations of parameters for the associated frequency of each mode. At the right column the period of oscillation of each mode assuming m = 2.5 × 10 − 22 ​ e ​ V 𝑚 2.5 superscript 10 22 𝑒 𝑉 m=2.5\times 10^{-22}eV and core radius r c = 1 subscript 𝑟 𝑐 1 r_{c}=1 kpc.

Mode	Detector $D_1$	Detector $D_2$	$\Delta$	T (Gyr)
$l=0$	${}_{}{}^{00}\nu _1^{}=0.0462$	${}_{}{}^{00}\nu _1^{}=0.0462$	0.00004	1.65
	${}_{}{}^{00}\nu _2^{}=0.0817$	${}_{}{}^{00}\nu _2^{}=0.0817$	0.0005	0.93
$l=1$		${}_{}{}^{10}\nu _0^z=0.0019$	$0.0002$	37.1
	${}_{}{}^{10}\nu _1^{\varrho }=0.0462$	${}_{}{}^{10}\nu _1^z=0.0462$	$0.0001$	1.65
	${}_{}{}^{10}\nu _2^{\varrho }=0.0722$		$0.0006$	1.05
		${}_{}{}^{10}\nu _2^z=0.0756$	$0.0002$	1
$l=2$	${}_{}{}^{20}\nu _1^{\varrho }=0.0462$	${}_{}{}^{20}\nu _1^z=0.0462$	$0.0001$	1.65
	${}_{}{}^{20}\nu _2^{\varrho }=0.0720$	${}_{}{}^{20}\nu _2^z=0.0721$	$0.0002$	1.05
	${}_{}{}^{20}\nu _3^{\varrho }=0.0909$	${}_{}{}^{20}\nu _3^z=0.0908$	$0.0003$	0.84
$l=3$		${}_{}{}^{30}\nu _0^z=0.0021$	$0.0001$	37.1
	${}_{}{}^{30}\nu _1^{\varrho }=0.04635$	${}_{}{}^{30}\nu _1^z=0.0463$	$0.0011$	1.65
	${}_{}{}^{30}\nu _2^{\varrho }=0.0720$	${}_{}{}^{30}\nu _2^z=0.0720$	$0.0008$	1.05
	${}_{}{}^{30}\nu _2^{\varrho }=0.0912$	${}_{}{}^{30}\nu _2^z=0.0912$	$0.0009$	0.84
$l=4$	${}_{}{}^{40}\nu _1^{\varrho }=0.0462$	${}_{}{}^{40}\nu _1^z=0.0462$	$0.0005$	1.65
	${}_{}{}^{40}\nu _2^{\varrho }=0.0719$	${}_{}{}^{40}\nu _2^z=0.0719$	$0.0002$	1.05

Equations

$$i\frac{\partial \Psi }{\partial t}$$

$${\nabla }^{2}V$$

$${\Psi }_0$$

$$\delta \Psi$$

$$a_{l0}$$

Full text

I Introduction

Ultralight spin-less dark matter is one of the dark matter candidates that shows properties that suit the current observational restrictions. Among the convenient properties of this particle are that it is ultralight, with mass of order $m\sim 10^{-22}eV$ Matos-Urena:2000 ; Hlozek:2014 ; Ostriker:2016 , which implies it has a de Broglie wave-length that can prevent structures from developing cusps at their centers, which may be a solution of the cusp-core problem (Chen:2016, ; Schive:2015, ; Du:2016, ; Robles:2013, ; Bernal:2016, ). At larger scales, there are now structure formation simulations showing that this candidate behaves just like CDM Schive:2014 ; Schive:2014hza . Another problem possibly addressed by this candidate is the abundance of substructure Matos-Urena:2000 ; Ostriker:2016 ; Marsh-Ferreira:2010 ; Hlozek:2014 ; Schive:2015 . An appealing property is that for bosons with such a small mass, the critical condensation temperature is high, of order $T_{c}\sim 1/m^{5/3}$ in the TeV scale, between the electroweak and quark epochs, at very early stages between $t\sim 10^{-32}$ and $10^{-12}$s after the Big Bang JuanMagana .

However, on local scales some interesting problems remain, such as the relaxation processes of structures that can involve the gravitational cooling process GuzmanUrena2003 ; GuzmanUrena2006 , dynamical friction Mocz2017 or appropriate damping terms ChavanisDamping . Nevertheless, structure formation simulations agree in that ultralight bosonic matter clumps into structures with universal profiles composed of a core, sometimes called solitonic profile, and a surrounding cloud of dark matter that suits a NFW profile Schive:2014 ; Schive:2014hza ; Schwabe:2016 ; Mocz2017 ; ChavanisCore .

Various oscillations of cored structures have been found to be interesting and potentially important, in particular the quasinormal modes of equilibrium configurations presented in GuzmanUrena2006 , which are associated to the oscillation of cores formed in structure formation simulations NewRef . The implications of these effects have been already translated into potentially observable effects on star clusters Marsh2018 . In a different context, the configuration resulting from head-on mergers has oscillations where the density changes orders of magnitude, an effect that eventually could affect the survival of luminous structures living inside such dynamical environment AvilopGuzman2018b . For instance, in GonzalezGuzman2016 it was shown that luminous matter, modeled as as an $N-$-body system within the gravitational potential of the condensate after a merger, might not survive such violent processes of relaxation and would be spelled out due to violent oscillations of the final configuration. In a less cataclysmic scenario, during the structure formation process, a core is constantly being perturbed by structures moving around it Mocz2017 . If some of the oscillations can reveal whether the model is viable or provide the fingerprints of cores, it is important to analyze how the latter respond to different types of perturbations. Oscillations in a way, may impose important restrictions or confirm the viability of the model, which is a strong a motivation to present a systematic analysis of the response of cores to perturbations in order to know the essential response modes that can eventually be identified with observations or predictions. Equally important is the construction of the oscillations spectrum of cores, which is important for their characterization. These are some reasons to carry out the present analysis, which focuses on the response of cores to various types of perturbations.

For this analysis we have two basic assumptions, the first one is that the dynamics of ultralight BEC dark matter structures is ruled by the GPP system of equations, which is the combination of the Gross-Pitaevskii equation for the condensate (Pitaevskii:2003, ), and Poisson’s equation sourced by the density of the condensate itself, whose solution provides the potential trap that contains the condensate bounded in a region. The second assumption is that cores are spherically symmetric equilibrium configurations of the GPP system GuzmanUrena2004 ; GuzmanUrena2006 as commonly assumed Schive:2014 ; Schive:2014hza ; Schwabe:2016 ; Mocz2017 .

In addition, one can take a perturbative approach since it is well known these equilibrium configurations are stable under a variety of spherical and axial perturbations, and they also show a late-time attractor behavior, which means that eventually fluctuations with an arbitrary initial profile would tend toward one of these equilibrium configurations GuzmanUrena2003 ; GuzmanUrena2006 ; BernalGuzman2006a . Therefore, the method used here to study the modes triggered by various perturbations, consists in the application of perturbations to equilibrium configurations, solve the full GPP system of equations in time and diagnose the evolution of physical quantities of the system. In order to be concise, in this paper we concentrate the analysis to axially symmetric perturbations only, and on the resulting oscillation modes that are excited.

The paper is organized as follows. In Section II we descrbe the perturbations used and the analysis in the frequency domain. In Section III the frequency modes are presented and in Section IV we draw some final comments and observations.

II Perturbations and numerical set up

The GPP system of equations written in variables absorbing the constants $\hbar$, $G$, the boson mass $m$ and the $s$-wave scattering length between bosons is

[TABLE]

where $\Psi=\Psi(t,{\bf x})$ is the parameter order describing the macroscopic behavior of the boson gas, $\rho=|\Psi|^{2}$ is the gas density, $\Lambda$ is the self-interaction coefficient, and the potential $V=V(t,{\bf x})$ is the gravitational potential sourced by the gas density itself, which in general can depend on time if $|\Psi|^{2}$ does. Here, as in most of the recent analysis of structure formation and interaction between structures, we assume $\Lambda=0$, which is the free field case commonly named the fuzzy dark matter regime Hu2000 .

II.1 Perturbing equilibrium configurations

These configurations are constructed assuming a harmonic time dependence of the wave function $\Psi=e^{-i\gamma t}\psi({\bf x})$ and spherical symmetry $\psi=\psi(r)$, where the label $e$ stands for equilibrium. In this case the equations above reduce to a stationary system because $|\Psi|^{2}$ is time-independent, which implies $V=V(r)$ is also-time independent and Schrödinger’s equation can be cast as a stationary equation with solution $\psi_{e}=\psi_{e}(r)$. The resulting equations are solved in the domain $r\in[0,r_{max}]$ as a Sturm-Liouville problem for the eigenvalue $\gamma$ as described in GuzmanUrena2004 , with boundary conditions of regularity at the origin for $\psi_{e}(r)$ and $V(r)$, isolation conditions at the outer boundary, namely a monopolar potential satisfying $-M/r_{max}$ and $\psi_{e}$, together with the requirement that its radial derivative should vanish within a given tolerance at $r_{max}$.

The solution defines a one-parameter family of solutions that relates the central value of the wave function $\psi_{e}(0)$ and the eigenfrequency $\gamma$. Due to the known scale transform that leaves the system (1) invariant, $\{t,x,\Psi\}\rightarrow\{\hat{t}/\lambda^{2},\hat{x}/\lambda,\lambda^{2}\hat{\Psi}\}$ for an arbitrary value of $\lambda$, if we solve for one central value $\psi_{e}(0)$ then the whole family for all other possible values of this central value will be automatically found. The workhorse configuration corresponds to the case $\psi(0)=1$, which in turn implies that the central density is also 1, has an eigenvalue $\gamma=0.69223$, mass $M=2.0622$ and $r_{95}=3.93$ GuzmanUrena2004 , all of which will be helpful later on. The scale invariance allows one to study the properties, evolution, and -for the purpose of this paper- the reaction to perturbations of only this workhorse configuration, and the results can be rescaled to all the other possible configurations. This is the reason why we will perturb this particular configuration in the reminder of the manuscript.

II.2 Perturbations

For our analysis, the perturbation applied to the wave function of the equilibrium configuration is a spherical shell with thickness $\sigma_{p}$, centered at a distance $r_{p}$ from the center of the configuration itself, launched with radial wave number $k_{r}$ and shape proportional to specific spherical harmonics. Since we only consider axially symmetric perturbations, we will only involve modes $(l,m)$ with $m=0$ and study the cases $l=0,1,2,3,4$. The perturbation profile is inspired by that in BernalGuzman2006b where it was used to show the stability of equilibrium configurations against non-spherical perturbations, and implemented also in AvilopGuzman2018b to understand an oscillation mode remaining after the head-on merger of two ultralight bosonic galactic cores. The perturbation profile has the following explicit form

[TABLE]

where $\psi_{e}$ is the wave function of the equilibrium configuration. In order to have a consistent gravitational potential, one solves Poisson equation sourced by the density $|\Psi_{0}|^{2}$ at initial time. Then we track the evolution of the perturbed system by numerically solving (1) and analyze its behavior.

We use perturbations proportional to each mode $l$ separately, in order to study how each mode affects the configuration. For this we ran a number of simulations with various values of the parameters in (2) with the aim of determining peaks in the frequency domain that are perturbation independent. Thus, for each value of $l$ we use shells centered at two different radii $r_{p}=2r_{95},~{}4r_{95}$, two different thickness of the shell $\sigma_{p}=1,~{}2$ and three different wave number values $k_{r}=0,1,2$. For each combination of these parameters, the amplitudes $A_{l0}$ were fine-tuned in such a way that the initial total mass of the perturbed configuration is $M=M_{e}+\delta M$, where $M_{e}=2.0622$ is the mass of the equilibrium configuration in standard code units (as in GuzmanUrena2004 ). Four values of $\delta M$ were used $\delta M=0.1\%,~{}0.5\%,~{}1\%$ and $1.5\%$ of $M_{e}$, where $\delta M>0$, meaning that the perturbation adds matter to the equilibrium configuration.

II.3 The code

In order to follow the evolution of the perturbed configuration, we solve the GPP system (1) with initial conditions given by (2). For this we use an upgraded version of the code introduced in BernalGuzman2006a ; BernalGuzman2006b ; this code uses cylindrical coordinates to solve the equations in the space-time domain $\varrho\in[0,\varrho_{max}]\times z\in[-z_{max},z_{max}]$ for $t\geq 0$ and two values of $\varrho_{max}=z_{max}$, specifically 20 and 40. These two domains suffice to contain multiples of the scale $r_{95}$ and more importantly, the use of two considerably different domain sizes with consistent results, as shown below, indicates the results are unaffected by boundary effects. The numerical implementation uses a finite differences method to discretize the system of equations on a two-dimensional uniformly discrete spatial domain, with resolution in a convergent regime. For the time integration we use the method of lines with a third order accurate integrator. More details about the numerical methods can be found in BernalGuzman2006a ; BernalGuzman2006b .

II.4 Boundary conditions

Boundary conditions are very important since we want to analyze oscillation modes because the wave function could be reflected from the boundaries and promote the formation of certain unphysical modes. In order to avoid this, for the wave function we impose a sponge consisting of an imaginary potential, which is non-zero near the boundary and has the profile defined in Refs. GuzmanUrena2004 ; BernalGuzman2006b . For the gravitational potential a multipolar expansion of the potential is used.

II.5 Detectors

It is very comfortable to analyze the evolution of the perturbed configuration by monitoring a scalar function of the density. For spherically symmetric scenarios $(l=0)$, its value at the coordinate origin is useful. Nevertheless, as we expect that the evolution can be more general in the case of non-spherical perturbations, we measure the value of the density at Detector 1 located at $D_{1}=(\varrho_{D},0)$ and at Detector 2 at $D_{2}=(0,z_{D})$ with the aim of measuring the oscillation modes along the two coordinate axes. Although we use detectors at various positions along the $\varrho$ and $z$ axes, there are no differences between the frequency modes, and we only present the results for the specific location with $\varrho_{D}=z_{D}=2\simeq r_{95}/2$.

II.6 Runs set up, output and analysis

In order to carry out a consistent analysis of the runs, the numerical parameters were standardized. The time domain for the runs is $t\in[0,3000]$ in code units, even though in some of the plots below a subset of this time domain is used. The production resolution is $\Delta\varrho=\Delta z=0.01\sim 0.025r_{95}$ and we use the same time step $\Delta t$ in all cases. The Courant-type factor $\Delta t/\Delta\varrho^{2}=0.2$, which implies $\Delta t=2\times 10^{-5}$. The output from the simulations are the time-series of the density at the coordinate origin and at $D_{1}$ and $D_{2}$. What is done next is to calculate the Fourier Transform (FT) of the time series with output in the frequency domain. In order for the FT to be accurate in the frequency domain, it sufficed the output of the time series to have a period of $1000\Delta t$. For a consistency check the FT was also calculated using a downsampling time-series of the output every $2000\Delta t$.

III Analysis and results

III.1 Spherical perturbations and consistency check

We apply a spherical perturbation using (2) for the case $l=0$, evolve the initial wave function $\Psi_{0}$, measure the density at detectors $D_{1}$, $D_{2}$ and calculate the FT of the resulting time series. These time-series measured by the detectors for a given time-window and their FT appear in Figure 1 for the particular case of $\sigma_{p}=1$, $\varrho_{max}=20$ and $\delta M=0.1\%$. As a consistency check, for this particular mode, we used an independent code that solves the GPP system in spherical symmetry, with only numerical perturbations as done in GuzmanUrena2006 , measured the density at the center as a function of time and calculated its FT. The purpose of this exercise is to verify that the peaks in the present analysis are reproduced with an independent code at least for spherical symmetry. The result is that there are two peaks indicating two dominating modes with frequencies at $\nu_{1}\sim 0.0462$ and $\nu_{2}\sim 0.0817$. The first of these frequencies was obtained using a perturbation theory analysis in GuzmanUrena2004 , which we will call fundamental mode, and the second one was numerically in BernalGuzman2006b within the context of the stability analysis.

III.2 Perturbations with $l=1$

In this case we follow the prescription (2) with $l=1$. The initial $\Psi_{0}$ is evolved using the three values of the wave number $k_{r}$, the four mass contributions of the perturbation $\delta M$, the two values of $\sigma_{p}$ and the two numerical domains described above. We produce the time-series corresponding to the density measured at detectors $D_{1}$ and $D_{2}$. Then we calculate the FT that we show in Fig. 2, for $\sigma_{p}=1$, $\rho_{max}=20$, the three values of the wave number $k_{r}$ and two out of the four mass contributions of the perturbation $\delta M=0.1\%,0.5\%$.

In the analysis, the frequencies we consider associated to the configuration and not to a particular feature of the perturbation, are those that show a peak for all the combinations of $\sigma_{p}$, $\varrho_{max}$, $k_{r}$ and $\delta M$, otherwise we consider it only a feature. Features are for instance the first peak measured by $D_{1}$ along the $\varrho-$axis with a frequency about $\nu\sim 0.004$, which is not a peak for $k_{r}=0$ and $k_{r}=2$. Another one is that measured by $D_{2}$ with frequency $\nu\sim 0.027$, because it is not a peak for any of the values of $k_{r}$ and the two values of $\delta M$ shown in Fig. 2.

With this in mind, we now turn to discuss the peaks we consider to be perturbation independent. The first important mode triggered by the $l=1$ perturbation, measured by detectors $D_{1}$ and $D_{2}$ with respective frequencies ${}^{10}\nu^{\varrho}_{1}={}^{10}\nu^{z}_{1}=0.0462$, corresponding to the fundamental mode. There are two other peaks measured by $D_{2}$: one with frequency ${}^{10}\nu^{z}_{2}\sim 0.0756$ that is triggered only by this perturbation and one is a peak with very low frequency ${}^{10}\nu^{\varrho}_{0}\sim 0.002$, that deserves some comments. This mode can be seen clearly from the time series of the density measured by detector $D_{2}$, as shown in Fig. 3 for the perturbation with $\sigma_{p}=1$, $\varrho_{max}=20$ and $\delta M=0.1\%$. Its period is of around 500 units of time, that corresponds to the first peak close to the axis in the right side of Fig. 2. This frequency mode is particularly interesting because it has been triggered only by $l=1$ and $l=3$ perturbations, and has been calculated in AvilopGuzman2018b with the aim of explaining long period modes excited during head-on core mergers.

III.3 Perturbations with modes $l=2$, $l=3$ and $l=4$

This time we use (2) with $l=2,3,4$. Again, the density is measured at $D_{1}$ and $D_{2}$ and the FT of the resulting time-series are shown in Fig. 4 for $l=2$, in Fig. 5 for $l=3$ and in Fig. 6 for $l=4$. These examples correspond to the particular case of $\sigma_{p}=1$, $\rho_{max}=20$, the three values of the wave number $k_{r}$ and two out of the four mass contributions of the perturbation $\delta M=0.1\%,0.5\%$. For all the wave numbers $\sigma_{p}$, $\varrho_{max}$, $k_{r}$ and $\delta M$ there are various frequencies that are triggered and are shown in the Table.

The frequencies in Table I were constructed as follows. For each value of $l$ a batch of runs was performed with all combinations of $k_{r},~{}\sigma_{p},~{}\varrho_{max},~{}r_{p}$ and $\delta M$, and for each run the FT was calculated. For each of these FTs the peak frequencies were found and averaged among the results in the batch. The number $\Delta$ is only the standard deviation calculated with each peak frequency within the batch.

IV Further discussion

IV.1 Spectroscopy and discussion

Collecting the information obtained, the dominant peak frequencies are summarized in Table I. The first interesting frequency is that of the spherical fundamental mode, which is excited by all types of perturbation with frequency $\sim 0.0462$. The perturbation with $l=3$ triggers this mode with an average value slightly off-set from the average value of the peak frequencies, but with contained within the uncertainty box.

The second most excited mode is that with frequency $\sim 0.072$, triggered by all perturbations except $l=0$. In the case $l=1$ it is measured only along the radial direction and along the two axes when applying perturbations with $l=2,3,4$. A specific property of the mode with frequency 0.072, excited by perturbations with $l=2,3,4$ is that it is the dominant peak when the perturbation corresponds to $l=2$ as seen in Fig. 4. This dominance happens not only for the specific parameters used to produce this figure, but also in all the runs carried out for $l=2$, with the various combinations of numerical and perturbation parameters. This is perhaps an important mode that can be triggered by the dynamics of matter, dark or luminous, with an $l=2$ distribution, and eventually usable to characterize cores.

A third and very interesting mode is the one triggered by $l=1$ and $l=3$ perturbations, which has low frequency $\sim 0.002$. This mode is measured along the $z$-axis. It has a frequency similar to that found in the density of configurations resulting from head-on mergers AvilopGuzman2018b . The distribution of matter in head-on mergers is consistent with the $l=1,m=0$ perturbation mode used here. This mode becomes important since head-on mergers result in a single blob whose density oscillates with amplitudes varying by 2 orders of magnitude. A catastrophic result of this scenario is that perhaps luminous matter could not survive such mergers as shown in GonzalezGuzman2016 , an effect that could serve to determine some predictions of this dark matter model or restrict in its viability.

A fourth frequency mode with frequency $\sim 0.091$ is found for $l=2$ and $l=3$ perturbations. This is the highest frequency peak described here. The reason for not exploring higher frequency peaks is that the resolution in the frequency domain becomes poor as frequency grows. Finally, there is a mode observed only for the spherical perturbation with frequency $\sim 0.0817$, discovered during the analysis of stability in BernalGuzman2006b .

IV.2 Translation to physical units

A useful method of unit translation is to look at the physical value of the equilibrium configuration eigenvalue $\gamma=0.69223$, for the workhorse configuration with central density equals to one GuzmanUrena2004 . This frequency corresponds to $\nu_{e}=\gamma/(2\pi)\simeq 0.11017$ in code units. In order to estimate the period associated with this frequency in physical units, one fixes the mass of the boson and the radius of the configuration as done in Mocz2017 , where the period corresponding to $\gamma$ or $\nu_{e}$ is given by $T_{e}=6.9\times 10^{8}\left(\frac{2.5\times 10^{-22}eV}{m}\right)\left(\frac{r_{c}}{1kpc}\right)^{2}$yr, where $r_{c}$ is the core radius that the density profile of an equilibrium configuration is fitted with, and $m$ is the boson mass. Using this formula with (for example) $r_{c}=1$kpc and $m=2.5\times 10^{-22}$eV, $T_{e}=.69$Gyr, with the translation of $T_{e}$ into physical units, it is now possible to translate the period of the oscillations triggered by the various perturbation modes, that we show in the right column of the Table. With these values of core radius and boson mass the time scales are of the order of Gyr. If smaller structures with $r_{c}=0.1$kpc are considered, the period of each mode will be two orders of magnitude smaller than those in Table I, which can have an effect on luminous matter on shorter time scales and thus impose restrictions on or make predictions about the model.

Finally a comment on the possible use of the spectrum of oscillations. The frequencies found here are the response to various types of perturbations. The implications of core oscillations are being investigated, for instance in relation to the possible formation of old star clusters near galactic cores. These analyses impose constraints to the mass of bosonic dark matter and represent state of the art studies Marsh2018 . So far those investigations are restricted to the diffusion regime, which is based on the consideration of a whole range of core oscillation frequencies. Now, most of the frequencies in our present analysis, except by the very low frequency ${}^{10}\nu_{0}^{z}$ mode, differ at most by only a factor of two, therefore the implications of analyses assuming a whole range of frequencies, would at the moment be unable to distinguish one mode from another mode in the spectrum presented here. However, with the results found in this paper, we expect that observable consequences of core oscillations can start being analyzed considering specific modes and more accurate constraints can be possibly found.

Acknowledgments

This research is supported by grants CIC-UMSNH-4.9, CONACyT 258726 and CONACyT 106466.