Compact objects in and beyond the Standard Model from non-perturbative vacuum scalarization

We consider a theory in which a real scalar field is Yukawa-coupled to a fermion and has a potential with two non-degenerate vacua. If the coupling is sufficiently strong, a collection of N fermions deforms the true vacuum state, creating energetically-favored false-vacuum pockets in which fermions are trapped. We embed this model within General Relativity and prove that it admits self-gravitating compact objects where the scalar field acquires a non-trivial profile due to non-perturbative effects. We discuss some applications of this general mechanism: i) neutron soliton stars in low-energy effective QCD, which naturally happen to have masses around 2 solar masses and radii around 10 km even without neutron interactions; ii) Higgs false-vacuum pockets in and beyond the Standard Model; iii) dark soliton stars in models with a dark sector. In the latter two examples, we find compelling solutions naturally describing centimeter-size compact objects with masses around 10^-6 solar masses, intriguingly in a range compatible with the OGLE+HSC microlensing anomaly. Besides these interesting examples, the mechanism of non-perturbative vacuum scalarization may play a role in various contexts in and beyond the Standard Model, providing a support mechanism for new compact objects that can form in the early universe, can collapse into primordial black holes through accretion past their maximum mass, and serve as dark matter candidates.

In this work, we consider a simple model-theory for NTSs in which a real scalar field h is coupled to a fermion through a Yukawa coupling.The scalar potential features the typical Mexican-hat shape, with two minima in h = ±v and a maximum in h = 0.However, the presence of the fermionic condensate S = ⟨ ψψ⟩ effectively modifies the shape of the scalar potential (see Fig. 2).In particular, if the Yukawa coupling f is sufficiently strong (or the scalar quartic λ sufficiently small), the point h = 0 becomes an actual (local) minimum of the theory.The latter condition is mathematically formulated by requiring that the parameter η = f /(2λ) 1/4 ≳ 1.In this regime, we show the existence of a bound NTS where the scalar field interpolates between h = 0 (false vacuum) and h = v (true vacuum) 1 , working in a fully relativistic approach.
The physical implication of the latter result is that, if η ≳ 1, a collection of N fermions is able to deform the 1 Without loss of generality, we are assuming that the asymptotic (scalar) ground state of the Universe is h = v.One can also consider a solitonic solution that connects h = 0 with h = −v, but then we would have a topological soliton (because the asymptotic value of the scalar field in the solution would be different from its value in the rest of the Universe).
FIG. 1: Sketch of the mechanism, showing the energy E of different configurations as a function of η, a fundamental parameter defined in Eq. ( 16).The bigger η the more the theory becomes strongly coupled.The orange balls represent massive fermions, whereas the yellow balls represent the false vacuum pockets of the scalar field in which the fermions are massless (and therefore depicted by smaller orange balls).(a) Standard ground state of the system.(b), (c) Whenever η ≳ 1, it is energetically convenient for the system to trap a fraction of fermions in false vacuum pockets.For η ≲ 1, equilibrium configurations describing bound false vacuum pockets do not exist.
h ≡ v ground state, giving rise to NTSs that describe false vacuum pockets, in which the N fermions are trapped (see Fig. 1).The latter are unable to escape (at least classically) because their energy is lower than N m f , the rest energy of N free fermions (m f is the fermion bare mass).Thus, the energy of an everywhere uniform configuration h ≈ v becomes higher than the energy of a configuration that shows false vacuum pockets here and arXiv:2401.06716v2[hep-ph] 7 May 2024 there.Since NTSs are intrinsically non-perturbative states (see Sec. III), the ground state of the system acquires a nontrivial scalar profile by means of non-perturbative effects.We dub this mechanism non-perturbative vacuum scalarization.
Investigations in a similar spirit were recenly carried out by Ref. [10] in the context of exotic neutron stars.
The primordial formation for NTSs has been investigated by many authors (see e.g.[11][12][13][14]), pointing out two main formation mechanisms, the solitosynthesis and the solitogenesis.In the former scenario, NTSs are formed by the fusion of N free fermions.This implies a non-zero cosmic asymmetry since it is necessary to accumulate a net number N of fermions over antifermions in a given region of space.In the latter scenario, a relic abundance of NTSs is produced through a first-or second-order cosmological phase transition.Recently, Ref. [14] reviewed these formation mechanisms in detail, showing that in certain cases NTSs can dominate the dark matter abundance.A further formation channel is provided from the Yukawa interaction present in Eq. ( 1), as long as it is enough long-range.In this case, fermions undergo clustering and structure formation even in a radiation-dominated era [15][16][17].
Phenomenological imprints of NTSs could be found through gravitational waves [18,19] and lensing observations.Besides being natural non-particle dark matter candidates, these solutions are also well-motivated examples of exotic compact objects [20] and might share similar phenomenology with ordinary compact objects such as black holes and neutron stars [21].

II. MODEL FOR RELATIVISTIC NTSS
Our starting point is the following theory, in which a (real) scalar boson h and a Dirac fermion ψ are minimally coupled to Einstein gravity2 , where the scalar potential is (see Fig. 2 for η = 0) Performing spontaneous symmetry breaking h → h + v in the Lagrangian (1) gives rise to a well-defined scalar  19)) normalized with respect to U 0 = λv 4 /16 as a function of h/v for various values of η in the Newtonian regime.In this regime, the zeros of the Fermi momentum are h/v = ±ω F (for |h|/v > ωF the fermions are removed from the theory).We fixed a representative value ωF = 0.5.We notice that when η ≳ 1 the effective potential develops a local minimum in h = 0.The shape of the potential is qualitatively the same also in the fully relativistic regime.mass term, However, for clarity of exposition, we prefer to work directly with (1), which has a simper analytical expression.The Yukawa interaction is controlled by the coupling f , giving an effective mass to the fermion, It is also useful to define which is the effective fermion mass when the scalar sits on the minimum v > 0.
The fermionic field has a U (1) global symmetry which ensures the conservation of the fermion number N (the Noether charge).
Our setup is similar to that of our previous work [22,23].An important difference is that here we relax the finetuning of the Yukawa coupling f imposed in [22,23], which allows us to explore connections with more realistic particle-physics content, as later discussed.
Fermions are treated through the Thomas-Fermi approximation [22][23][24], practically meaning that they enter Einstein's equations as a perfect fluid characterized by an energy-momentum tensor of the form where W is the energy density and P is the pressure of the fluid, while they also enter the scalar field equation through the scalar density S.These quantities are defined as follows where Notice that W = W (x µ ) through the spacetime dependence of k F and m eff (the same holds for P and S).The integrals in Eq. ( 8), ( 9), ( 10) can be computed analytically as shown, for example, in Ref. [22].
The fermion fluid is fully characterized once the Fermi momentum k F is given.Within the Thomas-Fermi approximation, it can be shown that where ω F is the Fermi energy at the origin (ρ = 0), which can be written in terms of the fermion central pressure P (ρ = 0) ≡ P c (see Ref. [22] for details).
In order to simplify the numerical integration, it is convenient to introduce the dimensionless quantities in terms of which the potential U and kinetic V = 1 2 e −2v(ρ) (∂ ρ h) 2 terms can be written as Moreover, we introduce the following dimensionless fermionic quantities Finally, the field equations (i.e. the Einstein-Klein-Gordon equations with the addition of the Fermi momentum equation) take the compact form [22,23] e where Ũ , Ṽ , P , W , and S depend on x, y, and r.The above equations contain two dimensionless combinations of the parameters3 , (16) the physical interpretation of which will be discussed in the next section.Here m p = G −1/2 is the Planck mass.
NTSs in the model (1) are static and spherically symmetric solutions to the above system of ordinary differential equations.The boundary conditions are the same as those used in Refs.[22,23], but in the present case4 y(r = 0) = ϵ, where ϵ is the initial displacement with respect to the false vacuum y = 0.The numerical procedure to build these solutions is the same as in Ref. [22], to which we refer for more details.
Considering more than one Dirac fermion or including Majorana fermions in Eq. ( 1) is a straightforward generalization.In the Thomas-Fermi approximation, this simply requires associating a Fermi momentum to each fermionic component and adding the relative energy density, pressure, and scalar density to the equations.

III. PROPERTIES AND INTERPRETATION
Before presenting the numerical results, it is useful to discuss some general properties of these solutions.
The NTSs arising from Eq. ( 1) are essentially fermion soliton stars with an effective positive cosmological constant inside, recently studied in [23].Therein, it was found that the energy of these solutions gets two contributions: one from the volume energy of the field and the other from a surface term.The relative importance of these two terms depends on the actual region of the parameter space.
In a field theory with spontaneous symmetry breaking, it is reasonable to expect v ≪ m p , which also guarantees that quantum gravity effects can be safely neglected.From Eq. ( 16), this implies Λ ≪ 1, which corresponds to the regime where the energy is dominated by the volume contribution (see Eq. ( 24) in Ref. [23]).In this regime both the mass and the radius scale as 1/Λ, while the parameter ωF scales as 1/η.Moreover, the binding energy Thus, for fixed Λ ≪ 1, the dimensionless binding energy is a linear decreasing function of η.Consequently, the bigger η the more bound the configurations5 , with η ≈ 1 marking the transition from bound to unbound states.As we will discuss in detail later, for η ≲ 1 the false vacuum pockets indeed disappear, and so does the nontrivial scalar profile of the ground state.Heuristically, this happens because the system has no benefit in investing some of its energy to build a false vacuum pocket to trap fermions that have a small mass gap between the two vacua of the scalar field.On the other hand, when η ≳ 1 the mass gap is big enough that it becomes energetically convenient to confine fermions in false vacuum pockets.For this reason, we will call the region confining regime 6 (in analogy to Refs.[22,23]).Moreover, the latter condition ensures that the point h = 0 in Eq. ( 2) becomes a local minimum.Indeed, the effective scalar potential arising from Eq. ( 1) is where U 0 = λv 4 /16.In order to analyze the behavior of the potential around the origin y = 0, we compute Thus, the point y = 0 will be a local minimum of the potential only if which in turn yields Using the fact that ωF ∼ 1/η, one gets again Eq. ( 18).In Fig. 2 we explicitly show the behavior of the effective potential for different values of η in the Newtonian regime, where the fermion scalar density S is analytically expressed as a function of y only.
It is important to stress that the NTSs discussed here are genuinely non-perturbative states.Indeed, in the limit in which one of the relevant couplings (λ or f ) is sent to zero, the energy diverges.When λ → 0, this happens because η → ∞.When f → 0, instead, the fermion mass vanishes (and so does the mass gap) and consequently the hydrostatic equilibrium implies h → v (see Sec. III A), as already anticipated from the previous heuristic argument.Therefore, NTSs arising from Eq. ( 1) actually become standard fermion stars [28] made of nearly massless fermions.In the zero-mass limit, the configuration does not have a finite radius and its energy diverges.In both cases, it is impossible to describe these solutions perturbatively around f ≈ 0 or λ ≈ 0.
Notice that if the two vacuum states are degenerate, whenever Λ ≪ 1 there is no need for a strongly coupled fermion since the confining regime is achieved as long as η ≳ Λ 1/2 (see Ref. [22]), so η can also be small.Indeed, in the latter case, the volume density of the scalar is zero, and the system needs much less energy to form a false vacuum pocket.
Finally, we highlight that our semi-classical approach, in which the scalar field is described as a classical solution (while quantum effects are taken into account for fermions), is valid whenever the scalar self-interactions are weak, i.e. λ ≲ (4π) 2 .

A. Numerical results
In Fig. 3 we show an example of a solution with the radial profiles for the metric, scalar field, and fermion pressure.
In Fig. 4 we present the mass-radius and compactnessmass diagrams for various values of η in the confining regime η ≳ 1.As anticipated, we observe that the massradius diagram is very weakly dependent on η.Moreover, the compactness becomes arbitrarily small along the tail of mass-radius diagram, giving rise to a Newtonian regime.
In Fig. 5 (left panel), we show the rescaled mass M Λ and rescaled radius RΛ, confirming the expected scaling in the Λ ≪ 1 limit.Overall, the mass-radius diagram is qualitatively similar to that of strange (quark) stars [29,30].
It is straightforward to identify a critical mass M c (and corresponding radius R c ) as the point of maximum mass in the M -R diagram.As confirmed in Fig. 5, for Λ ≲ 0.1 those quantities behaves according to As shown in Fig. 4, in the confining regime the dependence of the critical quantities on η is very weak.Thus, A and B are approximately constant, and numerical fits show that A ≈ 0.36 and B ≈ 1.35.We therefore obtain where q = (m h v) 1/2 , in terms of which the condition (18) can be written as We fixed Λ = 0.075.Varying Λ does not change the results qualitatively.There exists a turning point in the M -R diagrams at low masses, which cannot be seen from the figure, that proceeds towards the Newtonian limit of small M and large R, similarly to what already shown in [22].
Finally, the compactness of the critical configuration is GM c /R c ∼ A/B ∼ 0.27, slightly higher than that of a typical neutron star.
In Fig. 5 (right panel) we plot the mass as a function of ωF for various values of η at fixed Λ = 0.075, in the confining regime.In the low-mass region of these curves, where m h M/m 2 p ≈ 0, ωF → 1.In that limit ϵ → 1, the fluid energy density (together with its number density) vanishes and thus the false vacuum pocket is lost.As long as we depart from the m h M/m 2 p ≈ 0 branch, ϵ becomes a small number and the false vacuum pocket is recovered.From Fig. 5 (right panel) we also notice that, for any given η, ωF cannot be arbitrarily small.This means that configurations in which ϵ is small (i.e.truly false vacuum pockets) exist only above a minimum value ωmin F .We numerically verified that the scaling for ωmin Notice that this quantity is independent of Λ.Since the number density n in the core is proportional to ω 3 F (see Sec. II C in Ref, [23] for further details), we estimate The latter quantity can be interpreted as the minimum Noether charge per unit volume that allows for the existence of false vacuum pockets in the model (1).The configuration with the minimum possible Noether charge density has some peculiar properties.Calling M min and R min , respectively, its mass and radius, numerical analysis shows that Hence, in the Λ ≪ 1 limit this special configuration is pushed towards the tail of the mass-radius diagram, where the compactness G M min /R min ∼ Λ 2/3 → 0. We conclude that NTSs arising from Eq. ( 1) scan a range of masses and radii starting from the configuration in Eq. ( 28) deep in the tail of the mass-radius diagram, where relativistic gravitational effects are negligible, up to the critical configuration in Eq. ( 24), above which the NTS is expected to undergo gravitational collapse.In Fig. 6, we show the initial displacement ϵ and the binding energy of the configurations as a function of η.As expected from the discussion in Sec.III, the bigger is η the lower is the binding energy.In particular, we observe that the binding energy of the critical configurations, marked with a violet dot in Fig. 6, decreases linearly in η, as expected.Around η ≈ 1, we find configurations with positive binding energy, in agreement with the estimate in Eq. (17).
From the left panel of Fig. 6, we observe that the bigger is η, the lower is the displacement ϵ, whereas, for η ≈ 1, the displacement is ϵ ≈ 1, meaning that the scalar field is practically already on the true vacuum.In the latter case, as anticipated from the discussion in Sec.III, the effective fermion mass is already very near to m f and the picture of the configuration as a false vacuum pocket is lost.Therefore, we confirm that the ground state acquires a non-trivial scalar profile only if η is above a particular threshold.The numerical analysis done within the Thomas-Fermi approximation shows that already with η ≲ 2 we enter the deconfining regime (with the appearance of several turning points in both the mass and the radius, similarly to what is shown in Ref. [22]) and for η = 1 all the configurations with E B < 0 have h ≈ v.

B. Stability of NTSs
There are three main independent mechanisms through which a NTS may decay in time.First of all, if the binding energy is positive, a configuration with N fermions trapped inside is triggered to disperse into free N particles by quantum and/or thermal fluctuations [1].As discussed in Sec.III A, as long as Eq. ( 25) is satisfied, the configurations are bound.
Second, NTSs may be unstable under classical perturbations.A full investigation of the latter point goes beyond the scope of this work.However, in the confining regime, the shape of the mass-radius diagram is qualitatively similar to that of strange (quark) stars [29,30].Therefore, we expect the configurations below the critical mass to be at least radially stable, as it happens for ordinary stars [31].
Third, the configurations could be unstable under fission into smaller configurations.Calling E(N ) the energy of a NTS with N trapped particles inside, whenever fission is forbidden.The above equation can be recast in the following condition (also known as the Vakhitov-Kolokolov stability criterion [32,33]) For the configurations laying below the critical mass, we numerically checked that Eq. ( 30) holds.

IV. APPLICATIONS OF THE MODEL
In this section, we discuss explicit realizations of the model under investigation, embedded in various theories.

A. Neutron soliton stars
A simple application of non-perturbative vacuum scalarization is in the context of the linear sigma model [34], a low-energy effective theory of QCD that implements spontaneous chiral symmetry breaking.In this context, the fermion in Eq. ( 1) is just a nucleon, while the scalar field is a scalar meson called σ.Using the benchmark values m h ∼ 500 MeV [35], v ∼ f π ∼ 130 MeV in Eq. ( 28), the minimal configuration is estimated to be while the critical one in Eq. ( 24) is with a compactness of M c /R c ∼ 0.27, that is slightly larger than that of an ordinary neutron star.Moreover, being the mass of the nucleon m f ∼ 1 GeV, we get η ≈ 4, ensuring that the configurations lay in the confining regime.These estimates are in agreement with the ones given in Ref. [3].
Remarkably, the heaviest neutron stars discovered [36] have a mass compatible with the estimate in Eq. (32).In the absence of a scalarization mechanism, it is well known that a standard degenerate Fermi gas of neutrons can support compact stars only up to ∼ 0.7M ⊙ [37], and nuclear interactions are therefore needed to explain heavy neutron stars.In our model, instead, we can achieve the same result by simply coupling the degenerate Fermi gas to a scalar field.
We are neglecting the energy density of the pions assuming that their mean values on the ground state are zero.However, they are inevitably present in the linear sigma model and the inclusion of a pion condensate could lead to interesting phenomenology [38].Exploring these effects, however, would require a more involved analysis.
Our simple estimate suggests that the ground state of the heaviest neutron stars can in fact be scalarized and their core is made of a false vacuum pocket where the chiral symmetry is unbroken [10,39].Such NTSs may be formed in the early universe, through the mechanisms mentioned in the introduction, or during the merger of two neutron stars, if the energy involved is high enough (∼ f π ) to allow for the formation of false vacuum pockets and consequent trapping of nucleons.

B. Higgs false vacuum pockets
Non-perturbative vacuum scalarization may play a role also in the Higgs sector of the Standard Model (SM).Indeed, neglecting gauge interactions, Eq. ( 1) naturally arises after electroweak symmetry breaking, as shown in Appendix A. In this scenario, h is the physical Higgs field and ψ a SM fermion.Using m h ∼ 125 GeV and v ∼ 246 GeV, the minimal configuration is found to be while the critical mass and radius are These configurations are compelling candidates for nonparticle dark matter and exotic compact objects.In-deed, it is intriguing to note that the critical mass range is naturally compatible with the OGLE+HSC observation of some microlensing events [40,41].However, nonperturbative vacuum scalarization occurs only if there is a strongly coupled fermion to the Higgs boson, that is, only if the condition in Eq. ( 18) is satisfied.For the Higgs parameters, Eq. ( 18) gives m f ≳ 175 GeV, which is very near the top quark mass m top .Nevertheless, as discussed in Sec.III A, numerical computations show that we need η ≳ 2 to form bound false vacuum pockets, meaning that m f should be at least ∼ 350 GeV.This suggests that nonperturbative vacuum scalarization does not occur in the SM 7 .In addiction, even if the quark top mass was above the threshold, the electric charge carried by the fermions would produce a repulsive Coulomb force that in general renders large configurations unstable [43,44].Indeed, the total electric charge is estimated to be Q ∼ (2/3 e) n R 3 .
For the minimal configuration this gives Q min ∼ 10 14 C, whereas for the critical one Q c ∼ 10 29 C. In both cases, the total electric charge is way above the maximum charge compatible with hydrostatic equilibrium [45], i.e. 10 −22 C for the minimal configuration, and 10 −7 C for the critical one.However, top quarks also carry color charges.This gives a possible caveat to the previous arguments and leaves open the possibility of forming Higgs false vacuum pockets using only SM fields.To support the latter hypothesis, it is enough to estimate the fermion number density n inside the object.Using Eq. ( 27), we obtain For the Higgs field parameters, n ∼ 7 × 10 47 cm −3 , which is orders of magnitude larger than the nuclear matter density, n nucl ∼ 10 38 cm −3 .Such high densities are expected to give rise to exotic states, in which interactions among fermions cannot be anymore neglected, such as a top/anti-top bound state [46] or a colored superconductor [47].These scenarios could allow for bulk matter neutrality with respect to both electric and color charges.A self-consistent computation where Eq. ( 7) is substituted by a new energy-momentum tensor, in which all these effects are taken into account, may lower the aforementioned threshold and allow for the existence of bound configurations.
Non-perturbative vacuum scalarization may naturally occur in extensions of the SM by a fourth generation of (heavy) fermions.However, the simplest models along this direction are ruled out [48].A viable possibility is adding only a fourth family of chiral leptons strongly coupled to the Higgs sector, without modifying the quark content.This would allow for evading bounds on the number of generations coming from the Higgs decay into gluons.
It is important to stress that, regardless of the particular beyond-SM context in which Eq. ( 1) is obtained, as long as h is the Higgs boson, there are no free parameters in the model, being the mass and radius of the configurations very weakly dependent on m f .

False vacuum pockets evaporation
Beyond the three standard mechanisms mentioned in Sec.III B, there could be other ways of destroying the object, depending on the specific embedding of the model.Here we mostly focus on the case of (compact) Higgs false-vacuum pockets, but the argument can be extended also to other models in which the scalar h has Higgs-like couplings.
Let us assume that the fermions interact with (asymptotically) lighter SM particles (other than the Higgs boson).Inside the pocket, being the fermions effectively massless, their lifetime is expected to be extremely long.However, on the boundary of the pocket, the fermions re-acquire their mass and decays into lighter SM particles become kinetically allowed.Under the assumption that all the decay products leave the pocket (effectively subtracting energy from the configuration), it is possible to give a rough estimate of the lifetime.The number of fermions in the boundary is h is the size of the region where the fermions acquire their mass and ω 3  F is the fermion number density.The number of fermions that leave the object per unit time is estimated as A naive guess is Γ(m f → something) ∼ m f .Moreover, each fermion decay subtracts from the object an energy ∼ m f .Therefore, the total rate of energy loss is The lifetime is estimated as We notice that in the limit N → ∞ the lifetime would be infinity.Now, Using the Higgs parameters, one gets Λ ∼ 10 −16 and The actual computation of the lifetime should be carried out using a background field method (see e.g.Ref. [49]) and is complicated by the fact that we are in a regime of strong gravitational field 8 .Moreover, there are several subtleties that can drastically change the latter estimate.First of all, our computation relies on the assumption that the fermions are a gas of weakly interacting particles.Conversely, if one adds sizeable interactions between fermionic particles, the lifetime may become much bigger, just like the lifetime of the neutron inside a nucleus could be bigger than the age of the Universe, whereas the lifetime of a free neutron is very short.Second, in a beyond-SM scenario, the decay into lighter SM particles can be forbidden if the fermions in the configuration possess extra symmetries (see e.g.[50]) or do not couple with other SM particles.
What can not be disentangled from the other SM content is the Higgs field.In principle, some Higgs quanta of the system may be converted into other SM particles, producing a net particle flux.Since the Higgs boson is in a non-linear wave configuration, computing the latter quantity is a non-trivial problem.
We give an estimate through the following argument.A NTS can be quantum mechanically described as a superposition of many coherent states (see e.g.[51]).Each soliton quantum represents an interacting state which has in principle nothing to do with the standard perturbative states of the theory.The soliton size R is set by the energy µ of these soliton quanta and the configuration is dominated by modes with momentum k ≲ µ.In particular, Since in our model R ∼ 1 m h Λ , it is natural to estimate µ ∼ m h Λ ∼ 10 −14 GeV = 10 −5 eV.Therefore, just by energy conservation, the soliton quanta are not able to produce asymptotically a SM (massive) particle.Even neutrino (asymptotically) production is forbidden since What we cannot exclude is asymptotically photon production.However, since the Higgs boson does not couple to photons at the tree level, the overall effect is suppressed by loop effects.A rough estimate for the number of photons produced N γ gives where α is the fine-structure constant and N µ is the number of soliton quanta which is estimated to be ∼ m 2 h v 2 R 3 /µ ∼ 1/Λ 4 .Since the energy of each photon is ∼ µ, the change in the total energy is roughly µ dN γ /dt, which in turn gives the time needed to destroy the object An important caveat of Eq. ( 44) is that we are estimating the matrix element Γ(µ → γγ) as if the soliton quanta were the standard perturbative Higgs quanta of mass µ.However, we have already highlighted that the soliton quanta are purely interacting states which have nothing to do with the standard perturbative states.Therefore, a proper computation is needed to reliably estimate Γ(µ → γγ).

Accretion-driven collapse in the early universe?
Even in those scenarios where the solitons evaporate (as in the case of the false-vacuum Higgs solitons), they might still play a role during the evolution of the Universe, for example, if they survive long enough to accrete a sizeable fraction of their mass.
Indeed, while for purely scalar solitons (such as boson stars) scalar accretion-driven collapse to a black hole is prevented by gravitational cooling [52][53][54], our solutions have a fermionic core which can undergo gravitational collapse upon accretion of ordinary matter.
In the simplest scenario, Higgs balls are produced in the radiation-dominated era with an initial mass M i .In principle, such objects are able to accrete significantly before evaporation and one can therefore estimate the time needed to reach the critical mass M c above which the object can collapse into a (primordial) black hole.Following Refs.[55,56], we estimate the mass M of the object at the time t as where t i is the time of soliton formation and t i ≳ GM i since M i must be smaller than the horizon mass t i /G.At t ≫ t i , the final mass asymptotically is and can be significantly higher than Assuming that the object decays in a time t d and requiring M (t d ) > M c for black-hole formation, we obtain where the last step is valid if t d ≫ t i .
The above condition might or might not hold in specific models.For the case of the Higgs ball, we assume M c ∼ 4 × 10 −6 M ⊙ , t d = t i + t decay , and formation at the electroweak phase transition (t i ≈ 2 × 10 −11 sec).Therefore, for η ≈ 1, these objects will collapse to black holes before evaporation only if which is slightly smaller than the maximum mass and than the horizon mass at the time of formation.

C. Dark soliton stars
The simplest beyond-SM scenario where nonperturbative vacuum scalarization is realized is in the context of the dark sector paradigm [57], where we interpret Eq. ( 1) as embedded in the dark sector.The SM Higgs field h SM interacts with the dark scalar h through the unavoidable scalar portal ∝ h 2 SM h 2 , h SM h 2 , h 2 SM h.We assume m h ≳ m h SM , in order to kinematically forbid the Higgs boson direct decay into the dark scalar, allowing for evading collider constraints on the listed couplings.
In this framework, the evaporation bounds discussed in the previous section can be easily evaded.Indeed, assuming that ψ is the lightest fermion in the dark sector, decays on the boundaries are forbidden.Moreover, remembering that NTSs arising from Eq. ( 1) are made of scalar quanta with characteristic energy µ ∼ m h Λ, as long as the soliton quanta cannot produce asymptotically SM Higgs particles.Combining the latter equation with the requirement m h ≳ m h SM , we finally get Since it is natural to assume Λ ≪ 1, as discussed in Sec.III, Eq (50) is easily satisfied even for a dark scalar much heavier than the SM Higgs boson.
If the scalar quartic λ is an O(1) number, the condition m h ≳ m h SM translates into q ≳ q SM ≈ 175 GeV.Therefore, we conclude that these dark soliton stars are expected to be stable and support compact objects of masses naturally in the subsolar regime, according to Eq. (24).In particular, for m h ∼ v ∼ 10 2 GeV, we have a stable configuration with parameters similar to Eq. (34).
The dynamical formation of NTSs arising in a similar framework has been recently studied in [9], through a firstorder cosmological phase transition.The configurations produced, dubbed by the authors Fermi balls, are the nonrelativistic limit of our solutions.If an initial distribution of Fermi balls is able to accrete and cool down, the final state will be a scalarized ground state (see (c) in Fig. 1), well-described by our solutions.

V. CONCLUSIONS
In this work, we outlined the non-perturbative vacuum scalarization as a general mechanism to support new solitonic objects that can form in the early universe and serve as dark matter candidates.For concreteness, we considered a theory with a real scalar field coupled to a fermion in the context of General Relativity and found solutions describing self-gravitating compact objects where the scalar field acquires a non-trivial profile due to non-perturbative effects.
Besides the specific examples discussed in this work, other scenarios where this non-perturbative vacuum scalarization mechanism can play a role are fourthgeneration models with extended Higgs sector [58], asymmetric dark matter models [59], mirror symmetries, minimal supersymmetric SM [60,61], Type II see-saw mechanism [62], grand unified theories, inflation and cosmological phase transitions [63,64].Moreover, the inclusion of gauge fields into the solutions could give rise the important effects (see e.g.Ref. [43,44,65]), which we plan to explore in future work.
Although NTSs are energetically favored and stable under (at least) radial perturbations, couplings to other fields or self-interactions might induce evaporation of the solution, the details of which depend on the specific embedding of the model in a given theory.Future work should assess the role of fermion interactions for the evaporation time scale in given models.Even in cases in which this time scale is short relative to typical astrophysical scales, in the early universe NTSs might have enough time to accrete past the maximum mass or merge with other objects, with the possibility of forming primordial black holes in both cases [17,[66][67][68][69][70][71].Interestingly, this scenario would produce primordial black holes with a mass fixed in terms of the maximum mass of the soliton, regardless of the formation epoch.As we have shown, certain realizations of our model naturally lead to compact objects in a mass range that is compatible with the OGLE+HSC anomaly for some microlensing events [40,41].
A further important extension concerns the dynamical formation of these solutions.Many different channels have already been proposed in the literature.For example, NTSs arising from Eq. ( 1) can be produced by a firstorder cosmological phase transition [9].Alternatively, the Yukawa interaction, if enough long-range, drives clustering and leads to the formation of compact NTSs [15][16][17].Further possible formation channels, worth exploring in future work, are the following.First of all, statistical fluctuations inevitably present even during a crossover, are in principle able to provide a large concentration of fermions [72] and can be, therefore, the dominant source of charge fluctuations which leads to NTS formation through solitosynthesis.Moreover, one could consider a configuration where there is a gas of N free (massive) fermions moving in the true scalar vacuum h ≡ v.If η ≳ 1 and for a sufficiently large perturbation of the scalar field with energy above a certain threshold, we expect to end up in the true (scalarized) ground state.
The idea that a system can scalarize non-linearly, i.e.only if perturbed above a certain threshold, has been already numerically studied for scalar perturbations around a Schwarzschild black hole in scalar-Gauss-Bonnet gravity theories [73].An interesting extension of our work would be performing similar simulations in our model, possibly within a cosmological scenario.

FIG. 2 :
FIG.2: Effective scalar potential (see Eq. (19)) normalized with respect to U 0 = λv 4 /16 as a function of h/v for various values of η in the Newtonian regime.In this regime, the zeros of the Fermi momentum are h/v = ±ω F (for |h|/v > ωF the fermions are removed from the theory).We fixed a representative value ωF = 0.5.We notice that when η ≳ 1 the effective potential develops a local minimum in h = 0.The shape of the potential is qualitatively the same also in the fully relativistic regime.

FIG. 3 : 30 FIG. 4 :
FIG. 3: Radial profiles of scalar field h, normalized with respect to v, metric functions u, v (left panel) and fermion pressure (right panel) for a typical configuration (Λ = 0.075, η = 4).Continuous lines represent numerical data, whereas dashed lines reconstruct the asymptotic behavior of the solution by fitting with the Schwarzschild spacetime.The mass and radius of this configuration are m h M/m 2 p ≈ 4.08 and m h R ≈ 14.96, the compactness is C ≈ 0.27, while the solution parameters are Pc = 8.8 × 10 −3 and log 10 ϵ = −19.00.The binding energy is m h E B /m 2 p ≈ −5.10.

4 FIG. 5 :
FIG.5: Left panel: Rescaled mass-radius diagrams for various Λ.We fix η = 4. Varying η does not produce appreciable changes as long as we stay in the confining regime η ≳ 1.As expected, there is an approximate scaling law as ∼ 1/Λ, which becomes more and more accurate as Λ → 0. Right panel: Mass as a function of ωF for various values of η.We fixed Λ = 0.075.Varying Λ does not change the results qualitatively.Notice the existence of a minimum value for ωF at any given η.

FIG. 6 :
FIG.6:The initial displacement (left panel) and binding energy (right panel) for various values of η.We fixed Λ = 0.075.Varying Λ does not change the results qualitatively.As expected, the bigger η, the lower the binding energy.The configurations corresponding to the critical ones are highlighted with a violet marker.