# Lattice QCD study of mixed systems of pions and kaons

###### Abstract

The different ground state energies of -pion and -kaon systems for are studied in lattice QCD. These energies are then used to extract the various two- and three- body interactions that occur in these systems. Particular attention is paid to additional thermal states present in the spectrum because of the finite temporal extent. These calculations are performed using one ensemble of 2+1 flavor anisotropic lattices with a spatial lattice spacing fm, an anisotropy factor , and a spatial volume . The quark masses used correspond to pion and kaon masses of and , respectively. The isospin and strangeness chemical potentials of these systems are found to be in the region where chiral perturbation theory and hadronic models predict a phase transition between a pion condensed phase and a kaon condensed phase.

^{†}

^{†}preprint: JLAB-THY-11-1332

## I Introduction

An important goal of nuclear physics is to calculate nuclear properties and processes from Quantum Chromodynamics (QCD). The only known way to perform such calculations in an ab initio manner is using the non-perturbative formulation of QCD on a space-time lattice (lattice QCD). Much progress has been made in recent years in studying few hadrons systems, see for example the recent review of Ref. Beane et al. (2010). Many channels of meson-meson, meson-baryon and baryon-baryon interactions have been studied using the finite volume behavior of two-particle energy levels Hamber et al. (1983); Lüscher (1986a, b, 1991). Ground state energies of three- Beane et al. (2009a) and four- Yamazaki et al. (2010) baryon systems have been computed for the first time in QCD and quenched QCD, respectively. To become a central part of nuclear physics, the successes of lattice QCD in the realm of few-hadron systems must be translated into the realm of many hadrons and the complexity frontier must be addressed. To this end, recent work has focused on numerical investigations in systems of up to twelve pions or kaons Beane et al. (2008a); Detmold et al. (2008a); Detmold et al. (2008b); Detmold (2008), and algorithms Detmold and Savage (2010) have been constructed to greatly extend these studies. An important outcome of these studies was the extraction of the two- and three- body interactions in pion and, separately, kaon systems. The shifts in the -meson energy levels at finite volume were calculated numerically and then compared with the analytical expectations derived in Refs. Beane et al. (2007); Tan (2008); Detmold and Savage (2008). These calculations also probed the relation between isospin(hypercharge) chemical potential and isospin(hypercharge) density, finding agreement with the predictions of chiral perturbation theory Son and Stephanov (2001a); Kogut and Toublan (2001).

In this work, we study more complicated multi-meson systems involving both pions and kaons. As in the single species case, the expected volume dependence of the energy of an -, - system has been computed Smigielski and Wasem (2009) and these calculations will enable us to extract the scattering lengths for the three two-body interactions, , , , and also to study the various zero momentum three-body interactions, parameterized as , , , (restricted to the maximal isospin in each case). Multi-meson systems such as these exhibit a rich phase structure Dashen and Manassah (1974); Son and Stephanov (2001a) and are phenomenologically relevant in a number of settings ranging from RHIC Klein et al. (2003) to the interiors of neutron stars Kaplan and Nelson (1986). At low temperatures, the ground-state of a finite density system of, for example, pions is expected to undergo a transition from a Bose-Einstein condensed phase Sawyer (1972); Scalapino (1972); Baym (1973); Dashen and Manassah (1974); Campbell et al. (1975a, b) to a BCS-type superconducting phase as the isospin chemical potential is increased, and other interesting phenomena such as -condensation may also occur Voskresensky (1997); Sannino (2003); Aharony et al. (2008). Similar effects are expected in kaon systems. The phase structure of a mixed system of pions and kaons is a less well investigated question. Studies in the NJL (Nambu–Jona-Lasinio) model, Barducci et al. (2005); He et al. (2005); Warringa et al. (2005), random matrix theory Arai and Yoshinaga (2009), and SU(3) chiral perturbation theory (PT) with non-zero isospin and hypercharge chemical potentials Kogut and Toublan (2001), suggest that pion and kaon condensed phases compete in interesting ways for different ranges of isospin and hypercharge chemical potentials.

In the following, we focus on numerical analysis of lattice QCD realizations of these complex multi-meson systems and on the methodology needed to extract the two- and three-body interactions. We calculate the 90 possible two-point correlation functions involving -s and -s for all . These results are then analyzed to extract the ground state energies of the appropriate quantum numbers, taking care to account for thermal excitations in the lattice volume. Finally, the resulting large set of energies is used to extract the seven underlying interactions discussed above. The extracted energies also allow us to study the chemical potential of these systems. This study is performed using a single ensemble of 2+1 flavor anisotropic gauge field configurations generated by the Hadron Spectrum Collaboration Lin et al. (2009) and so should be viewed as an exploratory study, testing methods necessary for analysis of complicated many body systems. Future work with different lattice spacings, volumes and quark masses will allow contact with experiment Pislak et al. (2001, 2003); Adeva et al. (2005); Batley et al. (2006); Adeva et al. (2000), in the case of two-body interactions, and with other lattice studies.

## Ii Multi-hadron correlation functions

In order to extract the energies of the multi meson systems, we will study the behavior of two point correlation functions

(1) |

where the pion and kaon interpolating operators are defined in terms of quark fields as and , respectively. This system corresponds to a system of ’s and ’s with total momentum, .

It is instructive to first consider the form of such correlation functions which can be generically expressed as

(2) |

for some interpolating operators, , with
commensurate quantum numbers.^{1}^{1}1The creation and annihilation
interpolating operators are not in general Hermitian conjugates of
one another. In a system with infinite (Euclidean) temporal
extent, application of the transfer matrix formalism shows that

(3) |

where the sum runs over all eigenstates with the quantum numbers of the operators under study. The are the energies of eigenstates and . At large Euclidean time, the correlator is dominated by the ground state energy of the system, . However, in a system with a finite temporal extent, , a natural choice of temporal boundary conditions for quark and gluon fields, anti-periodic and periodic respectively, will force to be periodic (assuming the operators under study are bosonic, as in Eq. (1)). The expected form of the correlation function can again be inferred:

(4) |

with . The contributions on the right-hand side of Eq. (4) can be separated into two distinct classes. The first class is defined by choosing to be the vacuum state () and correspond to the ground state and excited states of the system. These contributions are special as they persist in the (or zero temperature) limit. The remaining contributions exist solely due to the finite time extent and are known as thermal states Beane et al. (2008a); Detmold et al. (2008b); Beane et al. (2009b); Richards (2001); Gockeler et al. (2002); Prelovsek and Mohler (2009); Detmold et al. (2009). These states vanish as because of the first exponential under the sum in Eq. (4).

Thermal contributions can be illustrated with a concrete example. Consider which corresponds to a system with the quantum numbers of four pions. The correlation function is dominated by a sum of three contributions in the large limit:

(5) |

where , the ellipsis denotes contributions involving excited pion systems that we will ignore in the current discussion, and the ’s are constant with respect to time. The first term in Eq. (II) corresponds to states where all four ’s propagate forward in time () and four ’s propagate backward in time (). The second term represents three ’s propagating forward in time and one propagating backward in time () as well as three ’s propagating backward in time and one propagating forward in time (). Finally, the last term arises from two ’s propagating forward in time and two ’s propagating backward in time ().

The general form for the correlation function of -pions and -kaons can be straightforwardly worked out and is given by the following

(6) |

where , the ellipsis denotes excited state contributions, and the last term in Eq. (6) is only present when are even. It is apparent that the number of possible terms contributing to grows with and . Locality of the transfer matrix guarantees that the eigen-energies appearing in many places in different are identical. Consequently, multiple correlation functions (choices of and ) can be used to extract the common set of eigen-energies.

Thermal states contribute to most lattice calculations and need to be considered in precise analyses. They are particularly prevalent in the multi-hadron context as these systems easily factorize into multiple color singlet states that can propagate over long distances. In Refs. Beane et al. (2008a); Detmold et al. (2008a); Detmold et al. (2008b), the issue of thermal states was avoided by using Dirichlet boundary conditions. In principle, such boundary conditions introduce unknown contamination into correlation functions, however these works analyzed correlations significantly separated from the boundary to reduce such effects.

## Iii Details of Lattice Calculation

The lattice calculations presented in this work are based on a single ensemble of anisotropic gauge field configurations generated by the Hadron Spectrum Collaboration Lin et al. (2009). These configurations were generated using a stout smeared Morningstar and Peardon (2004) gauge action and a 2+1 flavor clover fermion action; further details can be found in Refs. Lin et al. (2009). The spatial and temporal lattice spacings are and respectively, where is the renormalized anisotropy. The lattices have a volume of , corresponding in physical units to , and the light and strange quark masses are such that the pion and kaon have masses of and . These configurations have been extensively studied by the HSC and NPLQCD collaborations. For this study, an ensemble of 400 configurations were chosen from a long stream of 10,000 trajectories and are well separated in Monte-Carlo time such that autocorrelations are reduced.

Our analysis makes use of light and strange quark propagators generated by the NPLQCD collaboration. For each configuration, randomly positioned APE Best et al. (1997); Güsken et al. (1989) smeared sources were used to generate approximately 75 propagators. The EigCG inverter Stathopoulos and Orginos (2010) was used to perform the multiple inversions efficiently. These propagators were then APE smeared at the sink and combined into the correlation functions of Eq. (1).

Naively, the number of contractions involved in the correlation functions of Eq. (1) is enormous, . To compute the requisite contractions, we extend the methods developed in Ref. Detmold et al. (2008a) to the mixed species case. Following the manipulations of Ref. Detmold et al. (2008a) it is straightforward to show that for arbitrary matrices, :

(7) |

with

(8) | |||||

(9) |

where the correspond to the source and sink smeared quark propagators of flavor . Additionally, one can generalize the results of the single species case to show that in the mixed species case,

(10) |

where

(11) |

By expanding the exponential to a particular order in and , and equating Eq. (III) and Eq. (10), one can identify the function .

The computational complexity of these two-species contractions is
significantly more than that of the single species cases studied
previously. The reader can find an example, the
correlator, in Eq. (32). As in the single
species cases Detmold et al. (2008a); Detmold
et al. (2008b), high precision
arithmetic is required to perform these contractions correctly and
this is implemented using the arprec and qd libraries
Bailey et al. (2002). Explicit calculations for all
correlation functions show these correlators vanish to the requisite
precision and are an effective check of the correctness of our code.
The computational cost of computing all contractions for
and all is approximately twenty minutes on a
single core, compared to a few seconds for the single species case.
While not available at the time that the current calculations began,
the recursive constructions developed in Ref. Detmold and Savage (2010)
provide an alternate, and computationally more efficient, way of
generating the contractions. The results of both methods agree to
arbitrary precision.^{2}^{2}2In future calculations that extend the
number of pions and kaons beyond twelve, these recursive techniques
will be critical as the explicit code generated for all the
requisite contractions would be unmanageably large.

## Iv Volume Dependence of Multi-hadron energies

As discussed above, the calculations of the correlators in Eq. (6) determine the energies of the mesonic systems. These can in turn be used to determine the interactions through the well-known results of Lüscher Lüscher (1986a, b) in the two body case and the results of Refs. Beane et al. (2007); Tan (2008); Detmold and Savage (2008); Smigielski and Wasem (2009) in the many meson case where a perturbative expansion in the inverse volume is performed.

### iv.1 Volume dependence of two particle energies

An extraction of the scattering lengths from the single species and mixed species two particle systems provides a baseline reference with which to compare the results of the multi-particle analysis. In the center-of-mass frame, we define the interaction momentum, , from the energy shift . is determined from lattice calculations and one first solves for the momentum Lüscher (1986b). Lüscher’s formula relates this momentum to the phase-shift, , of two particle scattering as

(12) |

where and

(13) |

is a regulated three dimensional zeta function. By expanding the left-hand side of Eq. (12) using the effective-range expansion, , the scattering length, , can be determined. In the absence of interactions, possesses poles at for . With interactions these poles are shifted. One then computes as:

(14) |

Using this formula with the associated one-body and two-body energies will yield precise determinations of the scattering lengths as the one- and two- body energies are the most cleanly determined and the above results only omit exponentially suppressed finite volume effects.

### iv.2 Volume dependence of multi-meson energies

In Ref. Smigielski and Wasem (2009), the energy shift of a system -pions and -kaons in a finite volume from the corresponding non-interacting system was calculated. The shift is given by:

(15) |

with and and

(16) |

(17) |

where the parameters are related to the scattering lengths and effective ranges through Detmold et al. (2008a); Smigielski and Wasem (2009):

(18) |

It is the parameters that will be determined in the current lattice calculations. The four volume dependent (but renormalization scale independent) quantities characterizing the momentum independent three-body interactions are defined by ():

(19) |

and

(20) |

with . The functions , , , and along with the coefficients are defined in Ref. Smigielski and Wasem (2009). The finite parts of and are scheme dependent quantities where changes in the value will be compensated by changes in ; the numerical values for the Minimal Subtraction (MS) scheme are given in Ref. Smigielski and Wasem (2009). However, the are scheme independent and these are the quantities that will be determined during a lattice calculation. Furthermore, the three-body interactions in the and cases depend on the mass ratio, . Finally, in the limit of or , this result simplifies to the previously determined -boson case while the limit with and all interactions set to be equal it simplifies to the -boson case.

## V Analysis Strategies

After calculating the above correlation functions on the ensemble of gauge configurations, we obtain lattice measurements of in Eq. (1). Provided that we can reliably deal with the presence of thermal and excited states, these allow us to extract the ground state energies of the -pion, -kaon systems. By matching these measurements onto the analytic expectations of Eq. (IV.2), the various two- and three-body interaction parameters can then be extracted. There are many ways in which such an analysis could proceed. One could attempt to directly fit the correlation functions in terms of overlap factors, and the scattering parameters and three-body interactions by inserting the explicit form of the energies from Eq. (IV.2). However, we choose to perform this as a two stage analysis, first extracting the various energies by fits to the correlation functions using the model functions of Eq. (6) and then subsequently fitting these energies to the analytic forms of Eq. (IV.2) to extract the interaction parameters. The primary advantage of this approach to performing the fits is increased stability, however, care must be taken to preserve the significant correlations between the different correlation functions which are computed on the same set of gauge configurations.

### v.1 Wavefunction Overlap

As shown in Eq. (6), the -factors for each term are in general distinct. This makes for a difficult analysis since there are large numbers of linear and non-linear fit parameters that must be determined. In the non-interacting system, these -factors are simply related. In the single-species case, the relation would be such that is common to each term. The factorial simply counts the number of ways that -pions propagate forward and propagate backward. Generalizing this to the multi-meson case leads to . In an interacting system, however, this relationship does not hold. Nevertheless, in this work, we study pion and kaon systems with relatively weak interactions, and therefore we employ the ansatz above for the wavefunction factors. We stress this is a crucial assumption, allowing us to subsequently fit all 90 correlators. As a check, analyses were also performed allowing the -factors to remain distinct for as many correlators as possible. Results from both methods were found to agree within uncertainty, but a full analysis using independent ’s proved unfeasible.

An alternative method that can be utilized for these systems is the method of Variable Projection (VARPRO) Kaufman (1975); Fleming (2004) which eliminates linear fit parameters by performing their minimization analytically. This approach is useful in the case of unrelated ’s from Eq. (4). Additionally, when many energies exist to be fitted, convergence to the minimum of parameter space is not always guaranteed, particularly if the minimum is flat. To aid convergence, Bayesian priors Lepage et al. (2002); Howson and Urbach (2005) can be implemented so the physically reasonable region of parameter space is immediately tested. In an earlier stage of our analysis, the VARPRO method (with priors) was used for correlators with lower total particle number. The results were found to agree within uncertainty with the results using the ansatz above.

As a final check, the fitted energies were substituted back into the form for the correlator and compared against the original data set. Consistency was achieved within uncertainty. Therefore, for the remainder of this work, we assume the simplified form for the wavefunction overlap.

### v.2 Thermal effects of Correlation Functions

A particular multi-meson correlator depends not only on its ground state energy in the large- limit, but also upon fewer body energies in the thermal states that were discussed earlier. To account for this, we began by fitting the one pion(kaon) correlator to determine () with a bootstrap method as will be explained below. Once these values were known, we turned to the two pion(kaon) correlator. This correlator depends on () and (). Instead of refitting both energies, we used previous bootstrapped values of () in the two pion(kaon) fit function. This is reasonable since the cleanest signal for the energy will come from its respective correlator. To preserve all covariances in the data, the same bootstrap sample of gauge field configurations was used for all correlators. This sequential procedure was used to build up the analysis and ascertain all of the multi-meson energies. The fitted results for the energies and their statistical and systematic errors can be found in Tables 1 and 2.

### v.3 Statistical Analysis

The correlation functions span the time interval and a suitable subinterval, where the fit formula is applicable, must be chosen. As a first step, the (periodic) correlation function data was reflected around its midpoint in order to improve statistics. Effective mass plots for all correlation functions were constructed, allowing initial estimates to be made of where (in time) contamination from excited states ceases. These plots revealed that excited state contributions are suppressed by . Therefore, in an effort to be as general as possible, fits were run on all time intervals, , such that and . Using intervals beyond these values showed no improvement in fit results.

In this work, we utilized the bootstrap method of statistical resampling to estimate uncertainties. To construct the bootstrap ensembles, we first averaged over correlators from all sources on each gauge configurations, for each . Let be a given correlation function for some fixed number of pions and kaons, and , calculated on a given gauge field configuration at time . Then we denote this set as , with . The data covariance matrix is defined from the full data set as:

(21) |

where .

Next, we randomly sampled configurations from the full data set and constructed a bootstrap sample, of elements, for each time . This sample was averaged to yield and this entire process repeated times where , resulting in: . For each of the , a function was determined as:

(22) |

where represents the fit model in Eq. (6). Each of the functions were minimized and results recorded.

In order to assess the reliability of the fits, we define the goodness-of-fit, , where is the number of degrees-of-freedom, for each . Once a best fit energy is selected on each bootstrap sample, we determine a weighted average for the ensemble,

(23) |

with weights , and a weighted standard deviation,

(24) |

Above, and .

To assign a systematic uncertainty for the extracted energies, we return to the fits performed to the correlators. This is a necessary component of the analysis since there exists arbitrariness in the choice of which and which were selected for the correlation function fits as explained above. The systematic error on the best fit energy is defined as

(25) |

where are the energies extracted from using the time intervals and .

### v.4 Two and Three Body Parameters

The second stage of our analysis now focuses on the extraction of the scattering lengths and three-body coefficients from the measured energies using Eqs. (IV.2) and (IV.2). These determinations rely on a second bootstrap analysis involving a resampling of the extracted energies. The bootstrapping procedure for a specific correlation function yielded energies, and these formed the bootstrap samples for the extraction of the two and three-body parameters.

Once the best fit multi-meson energies were known, a very similar procedure used for the analysis of the correlators was used to find the ’s and ’s. Since a bootstrap ensemble exists for every best fit energy value, we created an energy sample, such that . This sample carries an additional vector index that labels the energies within the vector. In the case of single species pion energies (the kaon case is identical), an energy vector initially composed of was used to fit to and . We included another energy and refitted the interaction parameters and repeated this until all the energies were exhausted. In the multi-species case, a base set of along with ten randomly selected multi-species energies was created and fits performed for all seven hadronic parameters. This first set thus made use of 34 different energies. This set was enlarged by one, the parameters were refitted, and the process repeated until all ninety energies were used. The energy covariance matrix used in these fits, is defined according to:

(26) |

such that , and the energy on each bootstrap is defined as:

(27) |

where is shorthand notation for the fit functions in Eqs. (IV.2)-(IV.2).

The systematic errors assigned to the ’s and ’s are more complicated than those of the energies. Given a particular energy set of energies that are used to make a determination of ’s and ’s there are different combinations of the intervals that must be fit in order to completely propagate the systematic uncertainties of the energies to those of the interaction parameters (it is because there is a for each best fit energy as well as its systematic counterparts corresponding to the shifted time interval in the forward and backward direction). Even in the single species case, when , there are already combinations. For the multi-species case, it is too costly to fit all these permutations. Rather, we only fit randomly chosen permutations and take the difference of the mean of this set from the best fit and as the systematic error. From fitting all permutations in the single species case, up to , it was seen the systematic error stabilized well before the total number of combinations was computed and we assume this is also the case for the two-species case.

## Vi Results

### vi.1 Energies

Using the methods discussed above, we extracted the energies of the mixed and pure species system, from all ninety correlators. The final extracted values are shown in Tables 1 and 2 below, along with their associated fit ranges. These energies are shown in a three-dimensional plot along with their respective uncertainties in Fig. 1.

The fits become progressively more difficult as the number of mesons grows because of the increasing thermal contamination. This is directly reflected in the quality of the fits decreasing for large meson number in both the pure species and mixed species case. Fits to example correlators are shown in Fig. 2.

### vi.2 Interactions

The extractions of interaction parameters from mixed meson energies were performed to yield the three scattering lengths and four three-body coefficients. This work builds upon the studies of Beane et al. (2008a); Detmold et al. (2008a); Detmold et al. (2008b) and presents the first measurements of , and since these parameters can only be measured within the framework of the mixed-meson system.

The most straightforward determination of the scattering lengths is given by using the eigenvalue relation from Eq. (12). Using this, we find,

(28) |

The three-body coefficients can only be determined within the framework of Eqs. (IV.2), (IV.2), and (IV.2). We also use this same analysis to provide a check on the above results. Given that our analysis provides multiple determinations of the interactions parameters for varying numbers of combinations of energies used in the fits, these must be combined in some way to obtain the final values. Since each separate extraction can be viewed as a somewhat independent measurement, the final value given is taken to be the mean from the set of all extractions. The final uncertainties on the extractions are combinations of statistical uncertainties, systematic uncertainties obtained from variation of the fitting windows as discussed in Sec. V.4, and a second systematic uncertainty determined from the standard deviation of the full set of extractions, combined in quadrature. The systematics are the largest source of uncertainties in the results. The individual extractions of the various parameters and the final extractions are shown in Figs. 3-7. The error bars shown combine the statistical, and systematic uncertainties as discussed in Sec. V in quadrature. The shaded regions with thin borders denote the final results and their uncertainties. For the mixed species extractions, the second shaded band with thick borders, denotes the range of uncertainty in the quoted values from the single species analysis. These are shown together so the reader can see the overlap region between both sets of results. The poorest behavior originates from where the mixed species results drift away from the pure species one. The final values of the interaction parameters for the single-species case are:

(29) |

whereas for the multi-species case we find:

(30) |

For the two-body parameters, the perturbative analysis (for both pure species and mixed species) is in agreement with the nonperturbative values. In the above results, only the final uncertainty, including all statistical and systematic contributions, is given as discussed above. A full list of our extractions can be found in Table 3.

In previous studies, NPLQCD Beane et al. (2006a), CP-PACS Yamazaki et al. (2004), and ETMC Feng et al. (2010) have measured the pion scattering lengths in the isospin, channel. The NPLQCD determination of the scattering length at , yielded whereas an analysis from multi-pion correlators Beane et al. (2008a); Detmold et al. (2008a) yielded an extraction of and (in these results, the first uncertainty is statistical and the other uncertainties are systematics as discussed in the original references. In the ETMC analysis of scattering, is measured at . Fig. 8 shows the dimensionless combination from the current work in comparison to the determinations by other groups at a similar pion mass. However, these results are at non-zero lattice spacing and correspond to different discretizations, so agreement is not necessary. The , scattering length was also determined by the NPLQCD collaboration in Ref. Beane et al. (2008b). An analysis of the two-point kaon correlator yielded a value of , again at MeV. Analysis of multi-kaon correlators Detmold et al. (2008b), led to and where the uncertainties are statistical and systematic respectively. The scattering length has been investigated in quenched Miao et al. (2004) and full QCD Beane et al. (2006b) and the unquenched determination at MeV is . Hence, it is clear the current results are generally consistent with other groups’ extractions. The mixed species three-body parameters are novel results and are found to be of natural size and positive.

### vi.3 Isospin and Hypercharge Chemical Potentials

As we have determined the dependence of the energy of the mixed meson systems on the number of pions and kaons, we can construct the isospin and hypercharge (or strangeness) chemical potentials using finite differences following Refs. Detmold et al. (2008b); Detmold et al. (2008a) where systems of pions and kaons were investigated separately. In Refs. Detmold et al. (2008b); Detmold et al. (2008a), remarkable agreement was found between the numerical results and the leading order (LO) PT prediction Son and Stephanov (2001b) of the relation between the isospin(hypercharge) density and chemical potential

(31) |

with . The situation here is more
complicated since there are finite differences acting in various
non-orthogonal directions; the differences between and
determine while linear combinations of
, and determine
^{3}^{3}3 with
and . One goal of this analysis is to see
where on the vs. phase diagram Kogut and Toublan (2001),
the states created in the lattice calculation lie.

In Ref. Kogut and Toublan (2001), leading order SU(3)
PT is used to predict three distinct phases for nonzero isospin
and hypercharge chemical potential. The first is the normal
phase where the ground state has a net particle number of
zero. The other two phases are the pion-condensed and kaon-condensed
phases. The transition between the kaon-condensed phase and the
pion-condensed phase is predicted to be a
first order phase transition, separated by the line ,
while the transition from the normal phase to either
condensed phase is expected to be of second order ^{4}^{4}4An
AdS/QCD based model Albrecht and Erlich (2010) finds these transitions to be of
first order. and are defined by the lines and . These predictions assume zero temperature and are likely softened
by the non-zero temperature at which the lattice calculation is performed Loewe and Villavicencio (2003).

In Fig. 9, both the lattice calculations of and the PT phase boundaries are shown (dashed lines). Data points corresponding to higher numbers of particle states are shown in a orange/reddish color, while lower numbers are given in a blueish/greenish color. Points with large uncertainties are excluded from this figure for clarity (the omitted data correspond to the highest particle numbers). It is striking that the calculated chemical potentials mostly lie near the first-order phase transition line predicted by PT. Further calculations with larger numbers of pions and kaons will be enlightening, but more complex probes of these systems may be needed to fully understand the states that have been produced.

## Vii Conclusions

In this work, we have numerically studied complex systems of mesons of two distinct flavors, like-charged pions and kaons, and used them to extract information about the two- and three- body interactions amongst pions and kaons. Where known, the interactions were found to be consistent with previous calculations, however, two mixed-species three-body interactions were determined for the first time. Additionally, the isospin and strangeness chemical potentials and phase structure of the system have been investigated, with the systems preferring to probe a region in the plane where PT predicts a first order phase transition.

A major aim of this work was to investigate technical issues that arise in the analysis of complex multi-hadron systems. Accounting for the thermal states that proliferate in such systems, which easily factorize into distinct color singlet states, proved challenging and future calculations should avoid this by using larger temporal extents. Additionally, a number of techniques to perform coupled fits to the correlators studied were investigated and found to be beneficial in the analysis.

In the future, calculations probing larger meson numbers will allow further investigations of the phase structure of these interesting QCD systems. To understand the structure of the condensed systems created in the current and future calculations, more complicated observables that access transport properties may be needed; investigations in this direction are under consideration.

###### Acknowledgements.

We would like to thank C. Aubin, J. Erlich, K. Orginos, A. Nicholson, M. Savage, Z. Shi, A. Walker-Loud, and J. Wasem for comments and many useful conversations, and the NPLQCD and Hadron Spectrum collaborations for the use of their propagators and gauge configurations, respectively. BS would like to acknowledge the NSC of R.O.C. WD is supported in part by US Department of Energy grants DE-AC05-06OR23177 (JSA) and DE-FG02-04ER41302 and by OJI grant DE-SC0001784 and Jeffress Memorial Trust grant J-968. Calculations made use of computational resources provided by the NSF Teragrid and DOE NERSC facility and local resources at the University of Washington and the College of William and Mary.## Appendix A Tables

(27,37) | ||

(27,37) | ||

(27,37) | ||

(27,37) | ||

(31,41) | ||

(31,41) | ||

(26,36) | ||

(31,41) | ||

(31,41) | ||

(31,41) | ||

(32,42) | ||

(25,35) | ||

(27,37) | ||

(27,37) | ||

(27,37) | ||

(27,37) | ||

(27,37) | ||

(27,37) | ||

(28,39) | ||