**3. Results**

To solve the eigenproblem of Equation (4) we employed the Scalable Library for Eigenvalue Problem Computations (SLEPc) [10], which is a software package for the solution of large sparse eigenproblems on parallel computers. The solver used for the eigenvalue problem is a Krylov-Schur solver whose implementation within SLEPc is suited for non-Hermitian problems. We compute about 150 eigenvalues of the non-hermitian Dirac operator, which are the closest (in modulo) to the origin.

As the simulations are carried out for physical pion masses, away from the chiral limit (i.e., *mud* = 0), we try to extrapolate the density *ρ*(*ν*) to *mud* + *i* · 0 rather than to zero neglecting, at first, possible corrections due to non-zero masses and temperatures. We evaluate *ρ*(*mud*) by using kernel density estimation (KDE), a non-parametric way to estimate the multivariate probability density function from the measured spectrum. Such a technique is implemented in the python library scikit-learn [11], which we employ for the analysis.

It can be observed, by inspecting the contour plots in Figure 2, how only for *μI* large enough, i.e., within the BEC phase, the spectrum is wide enough in the real direction to encompass the red dot in Figure 2 at *mud* resulting in *ρ*(*mud*) = 0. At *μI* < *<sup>m</sup>π*/2 the eigenvalues are, instead, clustered along the imaginary axis and *ρ*(*mud*) = 0. At the largest simulated *μI* values, there is a tendency *ρ*(*mud*) → 0 due to the drift of the eigenvalues away from the real axis. However it should be noted that the impact of cutoff effects for larger and larger *μI* values remains to be assessed by a systematic comparison with the results on finer lattices.

Quantitative results for the spectral density are shown in Figure 3. It is interesting to match the *μ<sup>I</sup>*and *T*- dependence of *ρ*(*mud*) with the location of the boundary of the BEC phase, as determined by the onset of the pion condensate Σ*π* and with the location of the deconfinement crossover within the BEC phase as hinted for by a specific value of the renormalized Polyakov loop, that is *P*ren. = 1. This is done both in Figures 3b and 4 by using results for Σ*π* and *P*ren. obtained in the same setup in [3].

What can be observed is that the signal for the extrapolated spectral density seems to become nonzero around *μ*BEC *I* (*T*), that is at the location of the BEC phase boundary for the considered temperature. However, results also show that the extrapolated spectral density drops to zero again at larger values of *μ<sup>I</sup>*. Notice that lattice artefacts are expected to suppress *ρ*(*mud*), just as they do with Σ*π* [6–8]. Disentangling the signal for the BCS-BEC crossover from discretization errors at large *μI* is therefore difficult and a more systematic study is certainly needed to draw realistic conclusions. As one can see from Figure 5, the extrapolated spectral density still shows significant volume effects as well as a dependence of the results on the pionic source *λ*.

**Figure 3.** *ρ*(*mud*) as a function of *μI*, as obtained at various temperatures on 16<sup>3</sup> × 6 (**a**) and 24<sup>3</sup> × 6 (**b**) lattices. For the latter case only two temperature values are displayed. The lower (upper) edge of the shaded areas is set by the *μI* value at which the pion condensate Σ*π* becomes nonzero (the renormalized Polyakov loop *P*r becomes 1) in the same setup.

**Figure 4.** *ρ*(*mud*) as a function of Σ*π* (**a**) and of *P*ren. (**b**). For the latter case the shaded areas correspond to the range of values *P*ren. takes within the green Bose-Einstein condensation (BEC) boundary in Figure 1b at the two considered temperatures.

**Figure 5.** (**a**) *ρ*(*mud*) as a function of *μI* for *Nt* = 6 and three different spatial volumes *Ns*. (**b**) The *λ* dependence of results for *T* = 155 MeV.

#### **4. Discussion and Conclusions**

The presented results clearly show that the extrapolated spectral density is sensitive to the BEC boundary. However, to be able to draw conclusions on whether, and in which *μI* range, there is sensitivity to the BEC-BCS crossover as well, a more systematic analysis is needed. Such an analysis will allow us to establish the expected quantitative connection between the measured density and the BCS gap. Larger volumes, finer lattice spacings, and a *λ* → 0 extrapolation must be considered, and this is ongoing work. Finer lattices, in particular, will help us in identifying lattice artefacts due to cutoff effects at large *μ<sup>I</sup>*. Moreover, given that the Banks-Casher relation that we intend to use as a prescription to connect the spectral density with the BCS gap is strictly valid only for *T* = 0 and in the |*μI*| <sup>Λ</sup>*QCD* limit, a generalization of this relation away from this limit is desired. In addition, we might have to consider larger isospin chemical potentials and smaller temperatures.

**Author Contributions:** All authors contributed equally and significantly in writing this paper. All authors have read and agreed to the published version of the manuscript.

**Funding:** The research has been funded by the DFG via the Emmy Noether Programme EN 1064/2-1.

**Conflicts of Interest:** The authors declare no conflict of interest.
