{'fraction_non_alphanumeric': 0.033712428919577576, 'fraction_numerical': 0.01380991064175467, 'mean_word_length': 4.524457741211668, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 4, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Structure and cooling of compact stars

We study the structure and evolution of neutron stars (NS) the interiors of which are modeled using microscopic approaches and constrained by the condition that the equation of state (EoS) of matter extrapolated to high densities should not contradict known observational data from compact stars and experimental data from heavy-ion collisions (HIC). We use modern cooling simulations to extract distributions of NS masses required to reproduce those of the yet sparse data in the Temperature-Age (TA) plane. By comparing the results with a mass distribution for young, nearby NSs used in population synthesis we can sharpen the NS cooling constraints.


Introduction
Among compact stars denoted here generically as NSs one can theoretically distinguish three main classes according to their composition: hadronic stars, quark stars (absolutely bare or with tiny crusts), and hybrid stars (HyS). This classification is based on the understanding of the EoS of matter under extreme conditions which is at present lively debated. Can astronomy of NSs together with experimental data from HICs, where the matter exceeds the nuclear saturation density allow us to identify the true "name" of a given observed compact star? To answer this "simple" question one needs to do a very complex job. We have to collect experimentally justified models of nuclear matter, extrapolate their properties to sufficiently high densities, make a compact object model based on that and return to observations for comparison of the results. Probably this cycle is only one loop of iteration like in any perturbation approach.
To be sure that one moves in the correct direction one needs to define constraints on our modeling from the investigation of high-density behavior of nuclear matter (NM) activated in such an experiments as planned when the new accelerator facility (FAIR) at GSI Darmstadt will be constructed. On the other hand the new observations in astronomy provide new limits for the mass and the mass-radius relationship, surface temperature etc. of compact stars. These data constitute stringent constraints on the equation of state of strongly interacting matter at high densities, see [1] and references therein. Of course experiments such as CBM (=compressed baryon matter) at FAIR or NICA at JINR which are dedicated for the investigation of the phase transition from hadronic matter to the quark-gluon plasma (QGP) in HIC have a natural interest for the theory of NSs and particularly for the question whether quark matter exists in the core of a NS.
Here in these lectures, we demonstrate that the present-day knowledge of hydrodynamical properties of dense matter allows to construct hybrid EsoS with a critical density of the deconfinement phase transition low enough to allow for extended quark matter cores and stiff enough to comply with the new mass measurements of compact stars [2]. As has been argued by Alford et al. [3] [4] there are several modern QCD-motivated quark matter EsoS which could provide enough stiffness of high-density matter to be not in conflict with the new mass and mass-radius constraints. When compared with hadronic EsoS it turns out that quark and hadron matter have rather similar hydrodynamical properties in the region of the deconfinement transition. This means that hybrid stars can masquerade as neutron stars once the parameters of a generic phenomenological quark matter EoS have been chosen appropriately [4,2].
However, we suggest a new tool for "unmasking" the composition of neutron star interiors which is based on the fact that the state of matter at high densities determines the statistics of both NS observables, the temperature-age (TA) data as well as the mass distribution.

Equation of state of neutron star matter 2.1 Hadronic EoS
There are numerous comparative studies of NM approaches for HIC and NS physics applications in which a representation of the NM EoS has been employed which is based on the nucleonic part of the binding energy per particle given in the form E(n, β) = E 0 (n) + β 2 E S (n), where β = 1 − 2x is the asymmetry parameter depending on the proton fraction x = n p /n with the total baryon density n = n n + n p . In Eq. (1) the function E 0 (n) is the binding energy in symmetric NM (SNM), and E S (n) is the (a)symmetry energy, i.e. the energy difference between pure neutron matter and SNM. Both contributions E 0 (n) and E S (n) are easily extracted from a given EoS for the cases β = 0 and β = 1, respectively. The parabolic interpolation has been widely used in the literature, see e.g. [5], and appears rather robust. From Eq. (1) all zero temperature EsoS of nonstrange NM can be derived by applying simple thermodynamic identities [6]. In particular, we obtain for the baryonic energy density ε(n) and pressure P (n) as well as the chemical potentials of neutron µ n (upper sign) and proton µ p (lower sign), respectively. The neutron star matter (NSM) has to fulfill the two essential conditions of β-equilibrium and charge neutrality where µ e and µ µ are the electron and muon chemical potentials, conjugate to the corresponding densities n e and n µ . In case of a first order phase transition a mixed phase could arise in some density interval, see [7]. In general, the local charge neutrality condition could be replaced by the global one. However, due to the charge screening for the phase transition to quark matter this density interval is essentially narrowed [8,9]. The effect of the mixed phase on the EoS is also minor. Due to Eq. (5) the chemical potentials for muons and electrons are equal, µ µ = µ e . The threshold of muons to appear in the system is determined by their mass. In NSM the baryonic and the leptonic parts are considered as an ideal mixture in β-equilibrium ε(n, β) = ε B (n, β) + ε e (n, β) + ε µ (n, β) , (7) P (n, β) = P B (n, β) + P e (n, β) + P µ (n, β) .
Under NS conditions one parameter such as the baryo-chemical potential µ b or the baryon density conjugate to it is sufficient for a complete description. In β-equilibrated neutron matter the chemical potential of the neutrons equals the baronic one, µ n = µ b . Applying Eq. (4) and Eq. (5) shows that the electron and muon chemical potential can be written as an explicit function of baryon density and asymmetry parameter, µ e (n, β) = 4βE S (n).
Both electrons and muons are described as a massive, relativistic ideal Fermi gas.
With the above relations only one degree of freedom, namely the baryon density, remains in charge neutral and β-equilibrated NSM at zero temperature. Within this description actual properties of NM depend on the behavior of E 0 (n) and E S (n) only. Both can be deduced easily from any EoS introduced in these lectures.
Here in the discussions we use the nuclear EsoS that originate from relativistic descriptions of NM, because in the range of densities relevant NSs and HICs relativistic effects are important.
A brief description of different approaches is given in the following. For more details, see [1] and references therein.
• Phenomenological models originate from [10] and are based on a relativistic mean-field (RMF) description of NM with nucleons and mesons as degrees of freedom.
The mesons couple minimally to the nucleons. The coupling strengths are adjusted to properties of NM or atomic nuclei. A scalar meson (σ) and a vector meson (ω) are treated as classical fields generating scalar and vector interactions. The isovector contribution is generally represented by a vector meson ρ. In order to improve the description of experimental data, a medium dependence of the effective interaction has to be incorporated into the model.
In many applications of the RMF model, non-linear (NL) self-interactions of the σ meson were introduced with considerable success and were later extended to other meson fields (see [11]). As an alternative to NL RMF models, appoaches with density-dependent nucleon-meson couplings were developed (see [12] and references within), where for a flexible description of the medium dependence several parameterizations were introduced. In the NLρ the isovector part of the interaction is described, as usual, only by a ρ meson. The set NLρδ also includes a scalar isovector meson δ.
The density dependent RMF models are also represented here by [13] (DD models). They are obtained from a fit to properties of finite nuclei (binding energies, radii, surface thickenesses, neutron skins and spin-orbit splittings). In the D 3 C model additional couplings of the isoscalar mesons to derivatives of the nucleon field are introduced.

• Dirac-Brueckner-Hartree-Fock (DBHF) approach
One of the microscopic approaches starting ab initio from a given free nucleon-nucleon interaction that is fitted to nucleon-nucleon scattering data and deuteron properties is the DBHF method [14]. In this approach the nucleon inside the medium is dressed by the self-energy Σ based on a T-matrix. The in-medium T-matrix which is obtained from the Bethe-Salpeter equation plays the role of an effective two-body interaction which contains all short-range and many-body correlations in the ladder approximation. It is possible to extract the nucleon self-energies from DBHF calculations which can be compared with the corresponding quantities in phenomenological RMF models, but this is not completely unambiguous as discussed in Ref. [15]. Here, we use recent results of (asymmetric) NM calculations in the DBHF approach with the relativistic Bonn A potential in the subtracted T-matrix representation.
• A bridging approach between fully microscopic and more phenomenological descriptions that adjusts the parameters of the RMF model to results extracted from microscopic approaches. As an example of this method, we use a nonlinear RMF model (KVR) (or slightly modified parameter set (KVOR)) with couplings and meson masses depending on the σ-meson field [16]. The parameters were adjusted to describe the SNM and NSM EoS of the Urbana-Argonne group [17] at densities below four times the saturation density. In these models not only the nucleon but also the meson masses decrease with increasing NM density. Being motivated by the Brown-Rho scaling assumption, see [18], and the equivalence theorem between different RMF schemes, these models use only one extra parameter compared to the standard NL RMF model (NL model).
The sets of parameters and references of used models are given in [1]. The variation in the NM parameters is directly reflected in the behavior of the energy per nucleon in SNM E 0 (n) and of the symmetry energy E S (n) at densities above saturation as shown in Fig. 1. The various models of this study predict considerably different values for E 0 (n) and E S (n) at high densities. Under the condition of β-equilibrium, however, the range of binding energy per nucleon E(n, β) shows a much smaller variation than expected from E 0 (n) and E S (n). This is shown in the right panel of Fig. 1.

Quark matter EoS
Since there are no solutions of QCD for quark matter in the nonperturbative domain close to the chiral/ deconfinement transition at zero temperature and finite density one needs an effective model for the quark matter EoS satisfying at least the symmetry requirements of QCD. In order to discuss whether deconfined quarks can exist in neutron stars a successful class of approaches is that of Nambu-Jona-Lasinio (NJL) type models having a Lagrangian with chiral symmetry which is dynamically broken in the nonperturbative vacuum. In our applications in astrophysics we use this microscopic approach, see [19] and references therein.
In present investigation we took a three-flavor chiral quark model with selfconsistently determined quark masses and pairing gaps [20] similar to the parallel developments in Refs. [21,22], but generalizing [20] by including an isoscalar vector meson current [2] similar to the Walecka model for nuclear matter. While a stronger scalar diquark coupling leads to a lowering of the phase transition density, an increase in the vector coupling makes the quark matter EoS stiffer.
The thermodynamics of the deconfined quark matter phase as described within a three-flavor quark model of Nambu-Jona-Lasinio (NJL) type can be investigated systematically by applying techniques of finite-temperature field theory [23]. The path-integral representation of the partition function is given by whereμ andm = diag f (m u , m d , m s ) are the diagonal chemical potential and current quark mass matrices. For a = 0, τ 0 = (2/3) 1/2 1 f , otherwise τ a and λ a are Gell-Mann matrices acting in flavor and color spaces, respectively. C = iγ 2 γ 0 is the charge conjugation operator andq = q † γ 0 . G S , η V , and η D determine the coupling strengths of the interactions.
The Lorentz three-current,q γq, vanishes in the static ground state of matter and is therefore omitted in the vector interaction term in (10), for further details see [2,20].
After bosonization using Hubbard-Stratonovich transformations, we obtain an exact transformation of the original partition function (10). In mean-field approximation, when the bosonic functional integrals are omitted and the collective fields are fixed at the extremum of the action we obtain Here, Ω l is the thermodynamic potential for electrons and muons, and Ω 0 is a divergent term that is subtracted in order to get zero pressure and energy density in vacuum (T = µ = 0). E a (p) are the quasiparticle dispersion relations. Here, ∆ AA are the diquark gaps.m * is the diagonal renormalized mass matrix and µ * the renormalized chemical potential matrix, . The gaps and the renormalized masses are determined by minimization of the mean-field thermodynamic potential (12), subject to charge neutrality constraints which depend on the application we consider. In the (approximately) isospin symmetric situation of a heavy-ion collision, the color charges are neutralized, while the electric charge in general is non-zero. For matter in β-equilibrium, also the electric charge is neutralized. Further details are given in Refs. [2,20,21,22].  We consider η D as a free parameter of the quark matter model, to be tuned with the present phenomenological constraints on the high-density EoS. Similarly, the relation between the coupling in the scalar and vector meson channels, η V , is considered as a free parameter of the model. The remaining degrees of freedom are fixed according to the NJL model parameterization in Table I of [24], where a fit to low-energy phenomenological results has been made.
For our investigations of hybrid star structure we use two hadronic matter EsoS, where DBHF is an ab-initio calculation for the Bonn-A nucleon-nucleon potential within the Dirac-Brueckner-Hartree-Fock approach [15], discussed in the context of compact star constraints in [1,25]. DD-F4 denotes a relativistic mean-field model of the EoS with density-dependent masses and coupling constants adjusted to mimick the behavior of the DBHF approach [13]. The transition to the quark matter phase is obtained by a Maxwell construction, where the critical chemical potential of the phase transition is obtained from the equality of hadronic and quark matter pressures. A discussion of the reliability of the Maxwell construction for the case of a set of conserved charges is given in [8]. The resulting phase diagram for the caseof symmetric NM, relevant for applications in HIC experiments as CBM @ FAIR, is shown in Fig. 2.

Constraints on symmetric matter from HIC
The flow data analysis of dense SNM probed in HICs [26] reveals a correlation to the stiffness of the EoS which can be formulated as a constraint to be fulfilled within the testing scheme introduced here.
The flow of matter in HICs is directed both forward and perpendicular (transverse) to the beam axis. At high densities spectator nucleons may shield the transversal flow into their direction and generate an inhomogeneous density and thus a pressure profile in the transversal plane. This effect is commonly referred to as elliptic flow and depends directly on the given EoS. An analysis of these nucleon flow data, which depends essentially only on the isospin independent part of the EoS, was carried out in a generic model in Ref. [26]. In particular it was determined for which range of parameters of the EoS the model is still compatible with the flow data. The region thus determined is shown in Fig. 3 as the dark shaded region. Ref. [26] then asserts that this region limits the range of accessible pressure values at a given density. For our purposes we extrapolated this region by an upper (UB) and lower (LB) boundary, enclosing the light shaded region in Fig. 3.
Its upper boundary is expected to be stable against temperature variations. The important fact is that the flow constraint probes essentially only the symmetric part of the binding energy function E 0 (n).
In the hybrid star two phases are separated with a jump in baryon density corresponding to derivative of the pressure with respect to the baryochemical potential. The same is shown also for isospin-symmetric matter in Fig. 3. A slight variation of the quark matter model parameters η D and η V results in considerable changes of the critical density for the phase transition. The appearent problem of a proper choice of these parameters we solve by applying the flow constraint shown in Fig.3. At first we fix the vector coupling by demanding that the high density behavior of the hybrid EoS should be as stiff as possible but still in accordance with the flow constraint. We obtain η V = 0.50 such that the problem of the violation of the flow constraint for the DBHF EoS at high densities is resolved by the phase transition to quark matter. The optimal choice for η D is thus between 1.02 and 1.03. In the following we will investigate the compatibility of the now defined hybrid star equation of state with CS constraints.

Compact star structure, Mass-Radius relation
The mass and structure of spherical, nonrotating stars, to which we limit ourselves in this paper, is calculated by solving the Tolman-Oppenheimer-Volkov(TOV)equation, which reads as where the gravitational mass m(r) inside a sphere of radius r is given by which includes the effects of the gravitational binding energy. Eq. (13) describes the gradient of the pressure P and implicitly the radial distribution of the energy density ε inside the star. The EoS, i.e., the relation between P and ε leaves as much unknowns as we need to solve the set of differential equations. We supplement our EsoS describing the NSs interior by an EoS for the crust. For that we use a simple BPS model [27]. Due to uncertainties with different crust models one may obtain slightly different mass-radius relations. The stellar radius R is defined by zero pressure at the stellar surface, P (R) = 0. The star's cumulative gravitational mass is then given by M = m(R). Each configuration is specified by the given central density value n(0) at radius r = 0.  [2]. Indicated are also the constraint on the mass from the pulsar J0751+1807 [28] and on the mass-radius relationship from the isolated neutron star RX J1856 [29] The dotted lines indicate the gravitational redshift, z = (1 − 2GM/R) −1/2 − 1, of photons emitted from the compact star surface. Present constraints on the mass-radius relation of CSs do not rule out hybrid stars.

Observational constraint on the Mass and Radius
On the plot of the Mass-Radius relationship in Fig. 4 we also show regions corresponding to different observed objects. These results constrain the model EsoS for both hadronic and hybrid star matter. These data are the following where j ≡ cJ/GM 2 ≃ 0.1 − 0.2 is the dimensionless stellar angular momentum. •

Mass-Radius relation from RX J1856:
The thermal radiation of nearby isolated NS RX J1856.5-3754 with blackbody spectrum temperature T ∞ = 57 eV [32] because of initially estimated distance (60 pc) has been considered as a self-bound strange quark star with radius R ∞ ≈ 8 km. For the new measurements (distance ∼ 117 pc), RX J1856 has to have a rather large radius of R ∼ 14 km (R ∞ = 16.8 km), when the mass is 1.4 M ⊙ [29].
The resulting lower bound in the mass radius plane is shown in Fig. 4. Particularly interesting is the constraint from RX J1856 and there are three ways to interpret this result: A) RX J1856 belongs to compact stars with typical masses M ∼ 1.4M ⊙ and would thus have to have a radius exceeding 14 km (see Fig. 4). None of the examined EsoS can meet this requirement.
B) RX J1856 has a typical radius of R ∼ 12 − 13 km, implying that the EoS has to be rather stiff at high density in order to allow for configurations with masses above ∼ 2 M ⊙ . In the present work this condition would be fulfilled for DBHF, DD and D 3 C. This M > 1.7 M ⊙ explanation implies that the object is very massive and it is not a typical NS which would have M < 1.5 M ⊙ , as follows from population synthesis models.
C) RX J1856 is an exotic object with a small mass ∼ 0.2 M ⊙ , which would be possible for all EsoS considered here. No such object has been observed yet, but some mechanisms for their formation and properties have been discussed in the literature [33].
An other constraint on the EoS of stellar matter could be the critical mass M DU corresponding to the threshold of the most effective cooling mechanism, so called Direct Urca (DU) processes n → p + e − +ν e , which is an ultimate cause of very fast cooling of stars when its mass is only slightly above the critical value M DU [34]. For each hadronic EoS the critical mass is shown in Fig. 5. It depends on the symmetry energy E S and occurs when the proton fraction exceeds 11% (without muons and about 14% with muons in the medium) [1].
To consider the DU as a strong constraint one needs to investigate the cooling problem together with the nucleon superfluidity which suppresses the cooling rates. Calculations show that indeed the values of the pairing gaps used in the literature are not sufficient to resolve the DU as a problem [35].

Cooling behavior of hybrid stars
The cooling behavior of compact stars belongs to the most complex phenomena in astrophysics. Therefore, the codes for its numerical simulation as developed   (13) and (14) for all EsoS introduced in Subsect. 2.1. Crosses denote the maximum mass configurations, filled dots mark the critical mass and central density values where the DU cooling process becomes possible. According to the DU constraint, it should not occur in "typical NSs" for which masses are expected from population synthesis [36] to lie in the lower grey horizontal band. The dark and light grey horizontal bands around 2.1 M ⊙ denote the 1σ and 2σ confidence levels, respectively, for the mass measurement of PSR J0751+1807 [28]. by a few groups contain inputs (cooling regulators) of rather different kind, see, e.g., [37,38,39,40,41]. Attempts to develop a Minimal Cooling Paradigm [38] by omitting important medium effects on cooling regulators [40,35] unfortunately result in inconsistencies and suffer therefore from the danger of being not reliable. To develop a paradigmatic cooling code as an open standard, however, is rather necessary to cross-check the present knowledge of the groups before more sophisticated mechanisms like anisotropies due to the magnetic field [37] or special processes in the NS crust or at the surface are taken into account. Therefore, it is still premature to attempt an identification of the NS interior from the cooling behavior.
In order to circumvent such a model dependence we employ a given cooling code developed in Refs. [39,42] and vary the matter properties such as EoS, superconductivity and star crust model such as to fulfill all constraints known up to now (mass, mass-radius, TA, brightness, etc.). Moreover, we try to use consistent inputs.
The main neutrino cooling processes in hadronic matter are the direct Urca (DU), the medium modified Urca (MMU) and the pair breaking and formation (PBF) whereas in quark matter the main processes are the quark direct Urca (QDU), quark modified Urca (QMU), quark bremsstrahlung (QB) and quark pair formation and breaking (QPFB) [43]. Also the electron bremsstrahlung (EB), and the massive gluon-photon decay (see [44]) are included.
The 1S 0 neutron and proton gaps in the hadronic shell are taken according to the calculations by [45] corresponding to the thick lines in Fig. 5 of Ref. [39]. However, the 3P 2 gap is suppressed by a factor 10 compared to the BCS model calculation of [45], consistent with arguments from a renormalization group treatment of nuclear pairing [46]. Without such a suppression of the 3P 2 gap the hadronic cooling scenario would not fulfill the TA constraint, see [35].
The possibilities of pion condensation and of other so called exotic processes are included in the calculations for purely hadronic stars but do not occur in the hybrid ones since the critical density for pion condensation exceeds that for deconfinement in our case [39]. While the hadronic DU process occurs in the DBHF model EoS for all neutron stars with masses above 1.27 M ⊙ , it is not present at all in the DD-F4 model, see the right panel of Fig. 6. We account for the specific heat and the heat conductivity of all existing particle species contributing with fractions determined by the β− equilibrium conditions. Additionally, in quark matter the massless and massive gluon-photon modes also contribute.
In the 2SC phase only the contributions of quarks forming Cooper pairs (say red and green) are suppressed via huge diquark gaps, while those of the remaining unpaired blue color lead to a so fast cooling that the hybrid cooling scenario becomes unfavorable [47]. Therefore, we assume the existence of a weak pairing channel such that in the dispersion relation of hitherto unpaired blue quarks a small residual gap can appear. We call this gap ∆ X and show that for a successful description of the cooling scenario ∆ X has to have a density dependence. We have studied the ansatz ∆ X = ∆ 0 exp [−α (µ/µ c − 1)], where µ is the quark chemical potential, µ c = 330 MeV. For the analyses of possible models we vary the values of α and ∆ 0 , given in the Table 1 of [48] and shown in the left panel of Fig. 6.
The physical origin of the X-gap is not yet identified. It could occur, e.g., due to quantum fluctuations of color neutral quark sextet complexes [49]. Such calculations have not yet been performed within the relativistic chiral quark models. The size of the small pairing gaps in possible residual single color/ single flavor channels [50] is typically in the interval 10 keV -1 MeV, see discussion in [51]. The specific example of the CSL phase is analyzed in more in detail in Refs. [52,53,54,55].

Cooling curves in the TA diagram
We consider the cooling evolution of young neutron stars with ages t ∼ 10 3 −10 6 yr which is governed by the emission of neutrinos from the interior for t < ∼ 10 5 yr and thermal photon emission for t > ∼ 10 5 yr. The internal temperature is of the order of T ∼ 1 keV. This is much smaller than the neutrino opacity temperature T opac ∼ 1 MeV as well as critical temperatures for superconductivity in nuclear (T c ∼ 1 MeV) or quark matter (T c ∼ 1 − 100 MeV). Therefore, the neutrinos are not trapped and the matter is in a superconducting state. In Fig. 6 we show the density dependence of the pairing gaps in nuclear matter [45,39] together with that of the hypothetical X-gap in quark matter [56,47,48]. The phase transition occurs at the critical density n c = 2.75 n 0 = 0.44 fm −3 .
In Fig. 7 we present TA-diagrams for two different hadronic models and in Fig. 8 for two hybrid star cooling models presented in Ref. [48]. In the cooling calculations presented the crust model, i.e. the T m − T s relationship between be temperatures of the inner crust and the surface [39], has been chosen as such to fulfill the TA test: each data point should be explained with cooling curve belonging to an admissible configuration.
The TA data points are taken from [38]. The hatched trapeze-like region represents the brightness constraint (BC) [57]. For each model nine cooling curves are shown for configurations with mass values corresponding to the binning of the population synthesis calculations explained in [48].
In [48] for these hybrid cooling scenario the logN-LogS distribution constraint has been considered and it has been suggested to use the marking of TA diagram with five grey values in order to encode the likelihood that stars in that mass interval can be found in the solar neighborhood, in accordance with the population synthesis scenario, see Fig. 8. The darkest grey value, for example, corresponds to the most populated mass interval 1.35 -1.45 M ⊙ predicted by the mass spectrum used in population synthesis.   On the TA-diagram we show also the other constraint given by the maximum brightness of CSs, as discussed by [57]. It is based on the fact that despite many observational efforts one has not observed very hot NSs (logT > 6.3 − 6.4 K) with ages of 10 3 -10 4.5 years. Since it would be very easy to find them -if they exist in the galaxy -one has to conclude that at least their fraction is very small. Therefore a realistic model should not predict CSs with typical masses at temperatures noticeable higher than the observed ones. The region of avoidance is the hatched trapezoidal region in Figures 7 and 8.
Altogether, the hybrid star cooling behavior obtained for our EoS fits all of the sketched constraints under the assumption of the existence of a 2SC phase with X-gaps.

Unmasking neutron star interiors using cooling simulations
Considering the interior of neutron stars as a "laboratory" where matter under conditions of extreme densities occurs [58] we hope that it is possible to identify the composition of neutron star interiors from their observable properties like, e.g., masses, radii, rotational and cooling evolution [59,60]. The question for the existence of quark matter inside a compact object has often been mismatched with that of the existence of strange stars, essentially made up of strange quark matter without a hadronic shell and thus very compact   with radii less than ∼ 10 km [61,62]. The measurements of high masses for objects like the pulsar PSR J0751+1807 [28] with M = 2.1 ± 0.2 M ⊙ or of the mass-radius relation from the RX J1856-3754 pointing to either large masses or large radii of R > 14 km [29] require a stiff equation of state and exclude standard models for hyperonic or quark matter interiors as well as mesonic condensates, see [63].
Therefore, all suggested signals which are based on large first order phase transition effects and changes in the mechanical roperties such as the timing behavior of pulsar spin-down [64], frequency clustering [65] or population clustering [66,67] of accreting NSs should not be applicable.
Here we discuss a new, sensitive tool for "unmasking" the composition of neutron stars which is based on their cooling behavior [68]. As the cooling regulators such as neutrino emissivities, heat conductivity and specific heat in quark matter might be qualitatively different from those in nuclear matter, due to the chiral transition and color superconductivity with some possibly sensible density dependence, the TA curves for hybrid stars could be significantly different from those of neutron stars.

Mass distribution from TA data
In order to reach the goal of unmasking the neutron star interior we use here a new method for the quantitative analysis of the cooling behavior consisting in the extraction of a NS mass distribution from the (yet sparse) TA data [68] and its comparison with the (most likely) mass distribution from population synthesis models of NS evolution in the galaxy [36]. This method can be described as follows. For a given cooling model defined by the EoS and the cooling regulators, to each configuration with a gravitational mass M corresponding cooling curve T (t; M ) can be determined. For each set of mass values M i , i = 0, . . . , N M which is defining the borders of N M mass bins as, for example, in the population synthesis one can mark a region with borders of pair of neighboring cooling curves T (t; M i ) and T (t; M i−1 ). This strip in the TA plane corresponds to the i th mass bin. As a measure for the number of cooling objects to be expected within this mass bin we chose where N cool denotes the total number of observed coolers used for the analysis and P j (T, t) is the probability density to find the j th object at the point (T, t) in the TA plane. In Ref. [68] the simplest ansatz has been made that P j (T, t) is constant in the rectangular region defined by the upper and lower limits of the confidence intervals corresponding to the temperature and age measurements, (T jl , T ju ) and (t jl , t ju ), respectively, . (17) Note that in the case when the exact age the object is known (e.g., for a historical supernova), the time-dependence of P j (T, t) degenerates to a δ-function and the t-integral in (16) can be immediately carried out, leaving us with a onedimensional probability measure. This method has been applied to the cooling models for hadronic and hybrid stars described in the previous section. The results for the extracted mass distributions are normalized to 100 objects, defining N (M ) = 100 N i /( NM i=0 N i ), and shown in Fig. 9.
As we see from Fig. 9, the results are very sensitive to the chosen cooling model. In the hadronic scenario the onset of the DU cooling mechanism drastically narrows the mass distribution around the critical mass for the DU onset, see Fig. 7. On the other hand the slow cooling model predicts more massive objects than could be justified from the independent population analysis.
When comparing the density dependence of the pairing gaps, given in the left panel of Fig. 6, with the extracted mass distributions for the corresponding hybrid models in the right panel of Fig. 9, the direct relationship between the superconductivity and the mass distribution becomes obvious.
The DU problem as it was previously discussed in the literature [39,1,69] was based on the intuitive understanding that the mass distribution can not be peaked at a critical mass value which accidentally is unique for all observed young objects. Our modification of the definition of the DU problem does not contradict that suggestion, but rather provides a quantitative measure which to some extent rehabilitates the validity of cooling scenarios including the DU process.
On the other hand, the EoS model should obey the mass constraints too. Therefore, using the models discussed in this work we demonstrate that the most preferable structure of the compact object is likely to be a hybrid star with properly defined color superconductivity of the quark matter state in the core.  Figure 9: NS mass spectra extracted from the distribution of cooling data for both hadronic EoS models (left panel) and for hybrid stars with X-gap models I-IV (right panel). For comparison, the mass disctribution of young, nearby NS from the population synthesis of Popov et al. [48] is shown at the bottom of the panels.

Conclusions
In these lectures we have discussed the recent developments towards a test scheme for the EoS of matter under the extreme conditions using observational data from compact stars and experimental results from heavy-ion collisions.
The application of this scheme to specific EsoS offers some interesting insights which reveal the discriminative power of their combined tests in a broad region of densities including the possibility of a deconfinement phase transition.
In addition to the astrophysical part of the tests the results for the elliptic flow [26] restrict the hardness of symmetric matter which is in opposition to the constraint from neutron star maximum mass and/or mass-radius measurements.
One of the main results concludes that no present phenomenological finding bears a strong argument against the presence of a QM core inside NSs [3,2].
The resulting phase diagram for stellar as well as symmetric matter including the transition from nuclear to quark matter and satisfying the constraints discussed here exhibits almost a crossover transition with a negligibly small coexistence region and a tiny density jump. Therefore, the differences between the mechanical properties of the hadronic and hybrid stars are very tiny and the problem of the cooling history of the compact stars becomes central for the discrimination between alternatives for the composition of their interiors.
We have presented mass distributions obtained by analysing cooling calculations for NSs and demonstrated that they provide a sensible measure for the composition of compact star matter, even if the mechanical properties of the compact objects are almost identical.