QCD at finite isospin chemical potential

We investigate the properties of QCD at finite isospin chemical potential at zero and non-zero temperatures. This theory is not affected by the sign problem and can be simulated using Monte-Carlo techniques. With increasing isospin chemical potential and temperatures below the deconfinement transition the system changes into a phase where charged pions condense, accompanied by an accumulation of low modes of the Dirac operator. The simulations are enabled by the introduction of a pionic source into the action, acting as an infrared regulator for the theory, and physical results are obtained by removing the regulator via an extrapolation. We present an update of our study concerning the associated phase diagram using 2+1 flavours of staggered fermions with physical quark masses and the comparison to Taylor expansion. We also present first results for our determination of the equation of state at finite isospin chemical potential and give an example for a cosmological application. The results can also be used to gain information about QCD at small baryon chemical potentials using reweighting with respect to the pionic source parameter and the chemical potential and we present first steps in this direction.


Introduction
The study of first principles QCD under extreme conditions is of vital importance to illuminate the properties of QCD dominated aspects of matter in the universe and for phenomenology to test and develop models which can be used to understand different aspects of matter surrounding us. Despite the advances and successes of simulating QCD on a lattice in the past decades most regions of the QCD parameter space are still mainly unexplored. Particularly challenging is the study of QCD at finite density, since it is affected by the well-known sign problem, hindering numerical simulations. In the sector of the light up (u) and down (d) quarks the finite density parameter space can be characterised by baryon, µ B = (µ u + µ d )/2, and isospin, µ I = (µ u − µ d )/2, chemical potentials, as conjugate parameters to the associated densities, n B and n I , in the grand canonical ensemble. The sign problem appears as soon the baryon chemical potential is non-vanishing, while QCD at finite isospin chemical but with µ B = 0 is sign-problem-free and thus permits simulations of lattice QCD.
While typical systems in nature involving strongly interacting matter, such as the early universe, compact stars or heavy-ion collisions, share both, non-vanishing baryon and isospin densities, the study of QCD at pure isospin chemical potential is important and interesting in its own right. The expected phase diagram based on the findings in chiral perturbation theory [1], depicted schematically Speaker, e-mail: brandt@th.physik.uni-frankfurt.de Speaker, e-mail: schmalzbauer@th.physik.uni-frankfurt.de  Figure 1. Left: Schematic picture of the QCD phase diagram at finite isospin chemical potential, based on the findings in chiral perturbation theory [1]. The different phases are explained in the text. The black line indicates the conjuctured pure gauge deconfinement transition extending into the space of finite µ I from the µ I → ∞ limit.
in fig. 1, is expected to develop phases similar to the ones in the phase diagram at finite baryon chemical potential. In particular, at zero temperature and small µ I the system shows the so-called Silver Blaze phenomenon [2], where the groundstate of the system is not affected by µ I . When the system crosses the threshold chemical potential, µ I = m π /2, charged pions can be created, leading to pion condensation [3,4]. The associated phase transition is expected to be of second order in the O(2) universality class. While unimportant for most of the situations with isospin asymmetric matter, pion condensation can potentially play an important role in the description of neutron stars and for nuclear physics. Pion condensation also goes hand-in-hand with a proliferation of low-modes of the Dirac operator, leading to numerical problems in the simulations. A similar accumulation of low modes is also expected at finite baryon chemical potentials above threshold. In the pion condensation phase simulations are only possible with the use of an infrared regulator in the form of a pionic source term with parameter λ which is introduced into the action, see [5,6]. Physical results are then obtained by extrapolating the results to λ → 0, providing the main challenge in the analysis step.
At small temperatures the above features of the phase diagram are expected to remain mainly unaffected up to a possible shift of the pion condensation phase boundary, see fig. 1. Around the chiral symmetry restoration/deconfinement transition temperature T c , the pion condensate is expected to 'melt' or 'evaporate'. Consequently, the pion condensation phase boundary can potentially be shifted to very large values of µ I when T > T c . It is then interesting to investigate the interplay between the phase boundary to the pion condensation phase and the chiral symmetry restoration crossover. For extremely large values of µ I a decoupling of the quark and gluon degrees of freedom and a firstorder deconfinement phase transition associated with the gluon sector of the theory is expected to take place [1]. Thus it has been concluded that there might be a second phase transition for large values of µ I , depicted by the solid black line in fig. 1, which would then end on a second order critical point [1].
In the past decade a number of groups have investigated the properties of QCD at finite isospin chemical potential on the lattice [5][6][7][8][9][10] and in a variety of other approaches (see for example refs. [11]). However, all of the studies so far have been done on coarse lattices using unphysically large pion masses and/or an unphysical flavour content. In [12] we have presented first results of our study in the setup with 2+1 flavours of stout-improved staggered fermions at physical quark masses. In particular, we presented a novel method for the λ → 0 extrapolation using the singular values of the massive Dirac operator and presented first results for our study of the phase diagram and the comparison to Taylor expansion from µ = 0 on N t = 6 lattices. In this proceedings article we update these results by showing new results for N t = 6 and first results from N t = 8 and 10 lattices. We also present first results for the measurement of the equation of state at finite µ I . The associated result for the pressure at T = 0 can be used to construct a hypothetical boson star made of pions. First results and prospects for such a construction will be presented in section 3.4. Finally, we show how our results can be used to obtain information about the phase diagram of QCD at finite baryon chemical potential using reweighting. A particularly interesting, and up to date mostly unanswered, question concerns the properties of QCD in the enlarged (µ B , µ I ) parameter space. Our µ B = 0 simulations are idealy suited to study this parameter space for small values of µ B . First results in this direction will be reported in section 4.

Simulation setup and λ-extrapolations
In this section we sketch the setup of the lattice simulation. In particular, we show how the pionic source term is included into the lattice action, define relevant observables and discuss the methods used for the λ → 0 extrapolations.

Lattice action
We consider lattice QCD with three quark flavors u, d, s at temperature T = (aN t ) −1 in a discretized volume V = (aN s ) 3 with lattice spacing a. The quark masses m u = m d = m ud and m s are chosen to be physical and tuned along the line of constant physics [13,14]. The partition function is given by a path integral over all possible gauge fields, where β is the inverse coupling, S G the Symanzik improved gauge action and are the Dirac operators in the light and the strange sectors, respectively. The quartic roots in eq. (1) originate from the rooting procedure to remove the unwanted tastes from the simulations which appear due to the use of the staggered Dirac operator / D(µ). The off-diagonal elements in M ud originate from the introduction of a pionic source term and break the residual U τ 3 (1) symmetry of the action at λ = 0 explicitly. η 5 = (−1) n x +n y +n z +n t is the staggered fermion equivalent to γ 5 . For λ ∈ R and m s > 0, Monte-Carlo techniques like RHMC can be applied directly because of following from

Observables and λ-extrapolations
To extract information about the phase diagram, we study the pion condensate, the light quark condensate and the isospin density, Plugging (1) into these definitions and making use of (3), they are explicitly given by The traces can be evaluated using stochastic estimators or in the basis of the singular values of the massive Dirac operator [12]. ψ ψ and π are subject to renormalisation and, following [12], we define the renormalised condensates as where we have introduced the pion mass m π = 135 MeV and the chiral limit of the pion decay constant f π = 86 MeV for the purpose of normalisation.
Measurements of the observables are done for the ensemble including the artificial pion source term proportional to λ. To obtain physical results it is thus necessary to extrapolate the results to λ = 0. This λ → 0 extrapolation is the crucial and most difficult step concerning the analysis of the data. As shown in [12] a naive extrapolation is cumbersome and can lead to large systematical uncertainties. An improvement program for the λ-extrapolations has been outlined in [12], based on the singular value representation of the operator traces mentioned above. The first step consists of a "valence quark improvement", where the observable is replaced by (an approximation of) its λ = 0 counterpart. The simplest way to do this would be to set λ = 0 on the right hand sides of (6). This, however, is not possible in practice, due to the accumulation of small eigenvalues at finite µ I , leading to problems with the inversions in the stochastic approximation of the traces. Furthermore, for the pion condensate the source term is explicitly needed to obtain a non-vanishing value. The singular value representation of the traces in eq. (6) provides an alternative for the computation of the traces and we have shown in [12] how it can be used to reformulate the pion condensate in terms of the density ρ(ξ) of singular values of the massive Dirac operator ξ (see [15] for the derivation in the massless case), π = π · ρ(0) .
A similar improvement can also be done for the other observables and the associated publication containing the details is in preparation. This valence improvement already removes most of the λdependence of the expectation value, as we will show in Sec. 4.1. The remaining λ-dependence can be further reduced, constituting the second step in the improvement procedure, by reweighting the resulting expectation value with the leading order expansion of the full reweighting factor in λ. The details are provided in [12]. The remaining λ-dependence is mostly flat and can be extrapolated to λ = 0 in a well controlled manner.

Thermodynamics at finite isospin chemical potential
In this section we will present our results regarding the phase diagram, the equation of state and compare our results to the ones from Taylor expansion around µ I = 0. In the following all of the results have already been extrapolated to λ = 0 using the machinery described above. We will also briefly discuss a possible cosmological application for the equation of state determined at T = 0.

Results for the phase diagram
We start by updating the results concerning the phase diagram at finite µ I shown in ref. [12]. In particular, we present new results for the crossover temperatures for µ I < m π /2 on the N t = 6, 8 and 10 lattices and the resulting phase diagram. We also discuss the location of the chiral symmetry restoration transition within the pion condensation region. The pseudocritical temperature of the crossover can be defined by the inflection point of the renormalised condensate from eq. (6), for instance. In terms of this definition the crossover temperature in the continuum at µ I = 0 is T c (0) = 155(3)(3) MeV [16]. In the present study we start with a slightly different definition and determine T c (µ I ) by requiring that Σψ ψ takes the value Σψ ψ T c ≈ −0.550 at T c (µ I ), see [12]. This value agrees with the result for Σψ ψ T c obtained in the continuum limit at µ I = 0 [16]. Note, that this definition is only valid as long as we are considering the crossover in the Silver Blaze region, µ I < m π /2, where the condensate at T = 0 is independent of µ I , but does not hold for µ I > m π /2.
In figure 2 we show the resulting phase diagram for N t = 6, 8 and 10. The results for the pion condensation phase are those from ref. [12] and have been extracted from the points where the system develops a non-zero pion condensate. As observed already in ref. [12] for N t = 6, the crossover temperature at µ I = 0 appears to lie somewhat below the temperature associated with the melting of the pion condensate at high chemical potentials. Both temperatures tend to increase slightly in the approach to the continuum, as expected for the crossover temperature at µ I = 0, which should approach T c (0) = 155(3)(3) MeV, but the qualitative picture remains unchanged. The chiral symmetry restoration transition, however, is a broad crossover, whereas pion condensation sets in via a true phase transition.The latter is supported by a finite size scaling study in ref. [12]. The crossover line shows a downwards trend for growing µ I , even though the associated curvature shows large fluctuations for a → 0. To answer the question whether the chiral symmetry restoration transition and the phase boundary of the pion condensation phase start to coincide for µ I > m π /2 we need to change the definition of the crossover temperature and determine T c via the inflection point of the condensate. While this is still work in progress (and will eventually also replace the definition for T c for µ I < m π /2), a look at the behaviour of the condensate for µ I > m π /2 can be suggestive for the results that one might expect. We show the condensate for different values of µ I > m π /2 on the N t = 6 lattice in fig. 3 in comparison to the condensate at µ I = 0. The different starting points of the condensate at T = 113 MeV are due to the fact that the T = 0 condensate changes its value in the pion condensation phase, i.e. for µ I > m π /2. Consequently the subtraction of the T = 0 and µ = 0 condensate -needed for the renormalization (7) -leads to curves that do not start at unity. The dashed vertical lines indicate the temperature of the pion condensation phase boundary for these values of µ I . There is a clearly visible trend that the inflection point of the condensate and the pion condensation phase boundary approach each other with increasing µ I , meaning that one would expect chiral symmetry restoration to set in at the boundary of the pion condensation phase.

Testing Taylor expansion
As outlined in the introduction, one of the main challenges for simulations in lattice QCD is the sign problem for non-zero baryon chemical potential µ B . For small values of µ B , the sign problem can be overcome either by reweighting (see section 4.2) or via the Taylor expansion method. In the latter, expectation values of observables are expanded around µ B = 0. The resulting expressions contain the derivatives of the observable with respect to µ B evaluated at µ B = 0, which can be computed numerically. The main problem of the method is the a priori unknown range of applicability for a fixed order of the expansion. A similar Taylor expansion can also be performed for non-zero isospin chemical potentials, so that our results can be used to explicitly check the range of applicability of the method.
As before (see ref. [12]) we will focus on the isospin density n I for which the Taylor expansion with respect to µ I is given by with the Taylor coefficients c 2 and c 4 (the expressions for c 2 and c 4 are provided in ref. [12]). To determine c 2 and c 4 we use a cubic spline interpolation of the Taylor expansion coefficients from ref. [17], where the same action has been used.
In figure 4 we show the results for the comparison for the N t = 6 lattice and two values of the temperature where we still reach within the pion condensation phase. Note that by construction, the Taylor expansion is expected to break down at the phase boundary. Indeed this effect is clearly visible for T = 124 MeV, while for T = 155 MeV the disagreement is less obvious. This may also be due to the fact that the data remains in the vicinity of the phase boundary due to the strong flattening visible in fig. 2. For both temperatures the data clearly follows the curves obtained from Taylor expansion to O(µ 3 I ), which for both temperatures becomes distinct from the curve obtained from O(µ I ) between µ I = 50 to 60 MeV with the present accuracy of the data (for T = 124 MeV the difference is only visible thanks to the improved λ-extrapolations, which have improved the accuracy for n I compared to the results presented in ref. [12]). In the left panel of figure 5 we show a similar comparison for T = 176 MeV, where we do not enter the pion condensation phase. From this plot we can see that the good agreement with Taylor expansion to O(µ 3 I ) extends all the way up to µ I = 200 MeV, at least. After that slight deviations seem to appear, signalling the importance of terms of O(µ 5 I ). To make this more quantifiable we show the lines of constant difference ∆ = | n I − n I Taylor NLO |, where the second expectation value is the one obtained from Taylor expansion to O(µ 3 I ), in the right panel of figure 5. The plot shows that, as expected, the data disagrees with the Taylor expansion as soon as one enters the pion condensation phase. Above the phase boundary to the pion condensation phase the good agreement of the data with the Taylor expansion at a fixed order (here O(µ 3 I )) extends to even larger values of µ I for larger values of T , consistent with the notion that Taylor expansion is actually an expansion in µ I /T . Apart from the results at N t = 6 shown here, results for N t = 8 are also available and show no significant deviation from the findings above. A more detailed study of the approach to the continuum is postponed to a forthcoming publication.

The equation of state
One of the main ingredients for theoretical studies of phenomena and objects in cosmology and nuclear physics is the QCD equation of state. It is used for the hydrodynamic modelling of heavy-ion Figure 5. Left: Same as in figure 4 for T = 176 MeV. Note the difference in scale and that there is no vertical line since we are not entering the pion condensation phase for this temperature. Right: Contour plot of the difference ∆ between the simulation results and Taylor expansion (see text). To generate the data for the plot, n I has been interpolated using a cubic spline (see also section 3.3).
collisions and for the construction of neutron stars, for instance, just to name two of its many applications. For these physical situations the main contribution to the equation of state comes from finite baryon chemical potential. Nonetheless, the equation of state is also affected by the presence of a finite isospin chemical potential, so that the associated effects should also be included for a complete description of systems with isospin asymmetry.
Here we focus on the equation of state at pure isospin chemical potential. The two main quantities we consider, and from which all of the other quantities can be computed, are the pressure and the trace anomaly In this proceedings article we will show first results for the pressure obtained from the N t = 6 lattices, while the computation of the trace anomaly and other observables, as well as the continuum limit are left for future publications. The pressure can be rewritten as where we made use of the fact that n I = ∂p/(∂µ I ). The main task for the computation of p(T, µ I ) is thus the computation of ∆p(T, µ I ), while p(T, 0) is known from the interpolation provided in [14], for instance. To compute ∆p(T, µ I ) we evaluate the integral using a cubic spline interpolation of the data for n I , shown for some temperatures in the left panel of fig. 6. The results for ∆p(T, µ I ) are shown in the right panel of fig. 6. We can see that switching on µ I leads to a general rise of the pressure which is, not surprisingly, stronger for low temperatures, where one enters the pion condensation phase.

An application: Pion stars
The availability of a first-principles hadronic equation of state is of tremendous importance in astrophysics, because it allows to study compact stars under realistic conditions. The phenomenon of pion condensation at T = 0 and thus finite values of the isospin density for µ I > m π /2 could, at least in principle, lead to the formation of cold and self-bound stars, consisting of, say, positively charged pions 1 . These stars would correspond to a class of hypothetical objects called boson stars [19] 2 . They are expected to be stable against pion decay, due to the massless nature of the pions in the condensate.
In this section, we show how to construct such pion stars with the equation of state obtained from our lattice simulations 3 . Since we are considering a cold star we should use the equation of state obtained at T = 0. For this case measurements with the desirable setup of physical quark masses and improved actions used above are still in progress. Here we use the data for the equation of state obtained from the setup used for the tests in the next section instead. The mass-radius relation of a static, spherically symmetric and relativistic star can be obtained by solving the Tolman-Oppenheimer-Volkoff (TOV) equation [23] dp(r) dr with the gravitational constant G = 6.70861(31) · 10 −39 (GeV) −2 , shell radius r, pressure p(r), energy density (r) and mass For a known relation (p), eqs. (13) and (14) can be solved simultaneously by numerical integration, starting from some central pressure p(0) up to the edge of the star at radius r = R, indicated by a dµ I n I (µ I ) and (r) = −p(r) + µ I (r) n I (µ I (r)) .
Eq. (13) can then be rewritten as and the solutions may be labelled by the central chemical potential µ I (0). To obtain an interpolation for n I (µ I ) which smoothly goes to zero at µ I = m π /2, we use a cubic spline interpolation of our lattice data matched to chiral perturbation theory (see [6]) around µ I = m π /2. The mass-radius relation obtained by solving eq. (13) for different initial values µ I (0) is shown in fig. 7. The masses of the resulting stars range up to around three solar masses, similar to the expected masses of neutron stars. At the same time the pion stars have a diameter which is slightly larger than for neutron stars. Note, however, that these results have been obtained with an equation of state describing positively charged pions. This leaves us with a highly charged star, which is unlikely to be stable when the effect of the electromagnetic interactions is included. Neutrality can be reinstated by including further charged particles and we are currently studying the resulting system.

Reweighting to µ B 0
Up to now, our investigations were restricted to vanishing baryon chemical potential µ B . In most of the systems of interest, however, the finite value of µ B plays a crucial role. Since direct simulations at finite µ B are hindered by the sign problem, one can only use indirect methods such as Taylor expansion or reweighting to obtain information about the enlarged phase diagram in the µ B − µ I plane. These methods are, in general, restricted to small values of µ B (for reweighting this is due to the well-known overlap problem), so that the phenomenologically interesting region beyond the threshold of nucleon production cannot be reached. Apart from this region at large µ B the enlarged µ B − µ I parameter space offers other conceptually interesting regions, such as the region around the µ B = 0 axis for finite µ I , where it is interesting to investigate the behaviour of the boundary to the pion condensation phase at finite µ B , see fig. 8.
Since the simulations include a pion source, i.e. are performed at finite values of λ, it is also necessary to reweight the data in λ apart from the reweighting in the chemical potential. In this section we introduce the methods used for both types of reweighting and present first preliminary results. For testing purpose we switch to a cheaper setup which has already been used in ref. [6]. The results shown here are obtained on 8 4 lattices for T ≈ 0, neglecting the s-quark contribution, implementing unimproved staggered fermions and using the simple Wilson plaquette gauge action. The lattice spacing is a = 0.299(2) fm and we employ bare quark masses of am ud = 0.025, resulting in a pion mass of about 260 MeV. Pion condensation sets in at around aµ I = am π /2 ≈ 0.2.

Reweighting in λ
The basic idea behind reweighting is the following. The expectation value of an observable in a target ensemble, for us the pure isospin ensemble at λ = 0 with partition function is rewritten as a reweighted expectation value in an auxiliary ensemble. Here the auxiliary ensemble includes the pionic source λ, so that with Above we have used M = / D(µ I ) + m ud as a short notation. The computational cost for the determination of the reweighting factors is immense, as we have to compute all singular values of the massive Dirac operator. The cost can potentially be reduced by using the leading order expansion of the reweighting factor (cf. ref. [12]) Here the pion condensate can be computed using stochastic estimators to reduce the computational effort. Obviously, using only the leading order term as the reweighting factor is an approximation and one needs to ensure that the associated systematic effect is below the statistical uncertainty. One of the major prerequisites to achieve this goal is a strong correlation between the two types of reweighting factors. In the left panel of fig. 9 we show a scatter plot of the results for the normalised full and leading order reweighting factors. The correlation between R λ and R LO is clearly visible in the plot, but we observed that independent of λ, there is a slight tilt in the correlation with respect to the R λ = R LO line, which becomes more severe the bigger µ I . This behaviour can be explained by the different response of the two reweighting factors to fluctuations of eigenvalues. The details will be discussed in an upcoming paper. In the right panel of fig. 9 we show the effect of the leading order reweighting compared to the full reweighting for the example of the chiral condensate for different values of λ. Also included in the plot are the results for the condensate evaluated at finite value of λ and those without reweighting, but the operator evaluated at λ = 0, as explained in sec. 2.2. The plot indicates that the systematic effect associated with the use of the leading order reweighting compared to the full reweighting factor is indeed below the uncertainties. In fact, once we have improved the operator the remaining effect is small. In addition, we see that for the shown values of λ we do not observe any overlap problem with the λ = 0 ensemble, since the reweighting from different values of λ leads to results which agree within uncertainties. The green squares correspond to the improved λ-extrapolations, discussed in sec. 2.2, and we can see that the associated λ → 0 extrapolation is flat and well under control.

Reweighting in µ
To reach into the µ B − µ I parameter space we still need to include the reweighting step with respect to the chemical potentials of the light quarks µ u and µ d . The partition function of QCD with two mass degenerate staggered quarks at arbitrary chemical potentials is given by Note that the factorization of the quark matrices is only possible due to the absence of a pionic source. Once more we can rewrite expectation values in this target ensemble as reweighted expectation values in the pure isospin ensemble with λ = 0 as The reweighting factor in (23) is complex as soon as µ u −µ d , i.e. when there are baryonic contributions µ B 0 to the quark chemical potentials 4 The complex nature of R µ reflects that direct simulations at baryon chemical potentials are not possible with standard Monte-Carlo methods.
Combining (22) with the λ-reweighting in (18), it is possible to compute expectation values for two-flavor QCD at arbitrary chemical potentials in terms of the isospin simulations including a pionic source as The observables O are measured in the target ensemble (21). Using the definitions in (5) together with the baryon density, the analogue of the isospin density, the important observables are given by 4 We write µ I instead of µ I because we do not necessarily keep µ I constant in the reweighting process. Note that in the target ensemble the pion condensate always vanishes due to the absence of a pionic source in (21). This is in general the case for λ = 0 and it is the reason why the associated spontaneous symmetry breaking cannot be observed in simulations at finite volume.
Our plan to obtain valuable information about the phase diagram is to measure these observables for many different (µ u , µ d ) pairs in the µ B − µ I plane using the combined reweighting procedure from (25). In principle this involves recalculating R µ , eq. (23), many times for each value of µ u and µ d , which would consume a lot of computer time. To circumvent this, we use a determinant reduction scheme, first presented in [24], to express the (µ u , µ d )-dependence of the fermion determinants in a closed analytic formula, It is then sufficient to calculate the eigenvalues p i of the matrix P just once per configuration. A detailed explanation how to construct P can be found in ref [25]. Using this determinant reduction, the computation of the observables in eq. (27) reduces to computing the derivatives of ln det M with respect to m ud numerically and the ones involving the chemical potentials as which follows directly from (28). For the numerical derivatives ∂ f (m)/∂m = lim ∆m→0 [ f (m + ∆m) − f (m)]/∆m we found the most stable behaviour for ∆m ≈ 10 −3 m ud . We show the application of the above machinery to the isospin density in fig. 10. We simulated at different values of the isospin chemical potential µ I ∈ {0.0, 0.1, 0.15, 0.18, 0.2, 0.22, 0.25, 0.3, 0.4} at λ = 0.0025 and performed reweighting as in (25) to an equidistant grid of µ u and µ d values. In order to adress the overlap problem, we chose the auxiliary ensemble such that the error on the reweighted observable is minimal.
The reweighted isospin density is observed to become nonzero for µ B > m π /2 − µ I , i.e. outside the triangle around the origin in the left panel of Fig. 10. The reweighted chiral condensate also changes markedly beyond the hypotenuse of this triangle, see the right panel of the same figure (note that here we did not include the renormalization factors from eq. (7)). The fact that this region is connected to the µ B > m π /2 line at µ I = 0, where the sign problem is known to be severe, calls for a critical interpretation of this tendency. We are currently increasing the statistics to see whether this behaviour remains and employ different methods to estimate the overlap.
Since we cannot access the pion condensate, we determine the phase boundary of the pion condensation phase indirectly via the chiral condensate. Since the chiral condensate is expected to decrease at the boundary of the pion condensation phase, a behaviour also seen in the simulations for µ B = 0, we can follow the line of the pion condensation phase by searching for the region where this decrease is visible in the figure. The chiral condensate decreases strongly in the darker region on the right of the plot and one can clearly see a bend towards larger values of µ I with increasing µ B .

Conclusions
In this proceedings article we have presented the current status of our study of QCD at finite isospin chemical potential with improved staggered fermions at physical quark masses. The most crucial step in the analysis, the extrapolation of the prefactor of the pionic source term λ → 0, has been done using the improvement scheme introduced in ref. [12], which we have briefly sketched in sec. 2.2.
We have presented new results for the phase diagram on N t = 8 and 10 lattices, supplementing the results for N t = 6 presented already in ref. [12]. In the approach to the continuum the qualitative features of the phase diagram remain unchanged. We observe pion condensation at small temperatures starting at µ I = m π /2 up to a temperature of about 90% of the crossover temperature (T c (µ I = 0) ≈ 155 MeV in the continuum), where the phase boundary starts to shift towards larger values of µ I , until it flattens out at around 1.1T c (0). The crossover line shows a slight downward curvature for smaller values of µ I . For µ I > m π /2 we are currently working on the determination of T c via the inflection point of the condensate. By comparing the behaviour of the condensate to the location of the pion condensation phase boundary we obtain first evidence that the chiral symmetry restoration temperature is consistent with the pion condensation phase boundary starting from µ I 0.6m π . Concerning the comparison to Taylor expansion around µ I = 0, we have extended our range in chemical potential for temperatures above the pion condensation phase boundary and found the range of applicability of the expansion to O(µ 3 I ) to increase with increasing temperature, in agreement with the expectations. Due to the increase in precision from the improved λ-extrapolations, we have found that for T < T c Taylor expansion to O(µ I ) fails to describe the data at µ I between 50 to 60 MeV, while the data remains consistent with the expansion to O(µ 3 I ) up to the phase boundary. We have also presented first results for the pressure evaluated at finite isospin chemical potential on N t = 6 lattices, which is the first step for our measurements of the full equation of state. In sec. 3.4 we discussed an application of our measurements of the pressure in terms of the construction of gravitationally stable pion stars and showed first, preliminary results for the resulting mass-radius relation. The results presented so far are unrealistic in the sense that the star is highly charged and the next mandatory step is the inclusion of further charged particles to obtain a neutral star.
For an exploration of the phase diagram in the µ B − µ I plane, we used a two-step reweighting procedure in the pionic source and the quark chemical potentials. We have shown the applicability of the leading order approximation for the reweighting in the pionic source for our test systems, which helps to reduce the computational cost of the reweighting procedure tremendously. We provided first results for small baryon chemical potential and obtained first evidence that the pion condensation phase boundary shifts towards larger values of the isospin chemical potential for increasing baryon chemical potential.