Chiral soliton lattice in QCD-like theories

Recently, it has been shown that the ground state of quantum chromodynamics (QCD) in sufficiently strong magnetic fields and at moderate baryon number chemical potential carries a crystalline condensate of neutral pions: the chiral soliton lattice (CSL). While the result was obtained in a model-independent manner using effective field theory techniques, its realization from first principles using lattice Monte Carlo simulation is hampered by the infamous sign problem. Here we show that CSL, or a similar inhomogeneous phase, also appears in the phase diagram of a class of vector-like gauge theories that do not suffer from the sign problem even in the presence of a baryon chemical potential and external magnetic field. We also show that the onset of nonuniform order manifests itself already in the adjacent homogeneous Bose-Einstein-condensation phase through a characteristic roton-like minimum in the dispersion relation of the lowest-lying quasiparticle mode. Last but not least, our work gives a class of explicit counterexamples to the long-standing conjecture that positivity of the determinant of the Dirac operator (that is, absence of the sign problem) in a vector-like gauge theory precludes spontaneous breaking of translational invariance, and thus implies the absence of inhomogeneous phases in the phase diagram of the theory.


Introduction
The concept of spontaneous symmetry breaking has been one of the most fruitful paradigms in modern quantum physics, underlying a broad range of fascinating phenomena from superfluidity to ferromagnetism and the existence of light pseudoscalar mesons in the hadron spectrum. Starting from the assumption that the ground state of a given quantum system has lower symmetry than its equations of motion, one can take advantage of the powerful tools of effective field theory (EFT) to work out a model-independent description of the system's low-energy physics in terms of the ensuing Nambu-Goldstone (NG) bosons. The development of low-energy EFT for broken symmetries was pioneered half a century ago by Weinberg [2] and others [3,4], and its applications are still far from being exploited.
Much less is known about the conditions under which a continuous symmetry of a quantum system can, or cannot, be spontaneously broken. With rare exceptions where absence of spontaneous breaking can be proven [5][6][7][8][9][10][11], one likewise simply assumes, based on experience or intuition, that the ground state of a certain system maintains the full symmetry of its dynamics. In this regard, an intriguing proposal was made in ref. [12], linking the absence of spontaneous breaking of spatial translations in vector-like gauge theories (hereafter simply QCD-like theories) to the positivity of the determinant of their Dirac operator. The latter condition is required in order for the standard importance sampling lattice Monte Carlo techniques to be applicable; its violation is usually referred to as the sign problem. 1 If true, the conjecture of ref. [12] would thus make a remarkable addition to the list of existing no-go theorems for spontaneous symmetry breaking and at the same time, on the practical side, seriously hamper first-principle studies of inhomogeneous phases in dense quark matter.
The latter are not of mere academic interest. Namely, it is now widely believed that a substantial portion of the QCD phase diagram is occupied by inhomogeneous phases of various kinds, see ref. [13] for a review. Apart from its intrinsic interest, spontaneous formation of nonuniform order also has important consequences for the thermodynamics of dense nuclear/quark matter. However, a majority of the predictions of inhomogeneous phases in the QCD phase diagram is based on model calculations neglecting the effects of order parameter fluctuations, which may eventually render such phases unstable at any nonzero temperature [14][15][16]. A rare exception of a phase whose prediction is modelindependent and under theoretical control, is the chiral soliton lattice (CSL). This is a macroscopic manifestation of the chiral anomaly, whereby matter with nonzero baryon density can be realized without baryons, through a solitonic condensate of neutral pions. The CSL state requires for its existence either a sufficiently strong external magnetic field [1] or global rotation [17].
The goal of this paper is to show that the CSL phase, or a similar inhomogeneous phase, also appears in the phase diagram of an infinite class of QCD-like theories free of the sign problem. This has two immediate consequences. First, we disprove the above-mentioned conjecture that absence of sign problem implies absence of inhomogeneous phases in the phase diagram [12]. Second, our model-independent prediction of the CSL phase in theories amenable to lattice Monte Carlo simulations might provide a useful test bed for future numerical studies of dense quark matter. (To the best of our knowledge, the first attempt to study a crystalline phase of baryonic matter in lattice simulation was made in ref. [18], based on the (1 + 1)-dimensional Gross-Neveu model, followed by refs. [19,20].) The plan of the paper is as follows. In section 2, we introduce the class of (pseudo)real QCD-like theories and construct the corresponding low-energy EFT at the leading order of the derivative expansion. The anomalous Wess-Zumino (WZ) term, which is crucial for the existence of the CSL phase, is discussed to some detail in section 2.1. The ground state and excitation spectrum of the EFT is then studied in the subsequent sections, first in section 3 in the limit of vanishing quark mass (the chiral limit), and then in section 4 for massive quarks. A discussion of semi-numerical variational minimization of the Hamiltonian is relegated to appendix A.
The main result regarding the structure of the phase diagram was already reported in ref. [21]. In this sequel, we provide additional calculational details, in particular the derivation of the WZ term in section 2.1 and the variational minimization of the Hamiltonian in appendix A. The analysis of the excitation spectrum in sections 3 and 4 is completely new, and provides further supportive evidence that the phase diagram reported in ref. [21] is correct, as well as an alternative signature of the presence of an inhomogeneous phase in the phase diagram.

Effective field theory setup
In this paper, we consider a class of QCD-like theories where quarks transform in a real or pseudoreal representation of the gauge group. Examples of the class of real theories include theories with an arbitrary non-Abelian compact symmetry group and quarks in its adjoint representation, or theories based on the G 2 gauge group [22]. A prominent example of the class of pseudoreal theories is the so-called two-color QCD, that is a theory based on the SU(2) gauge group with quarks in the fundamental representation. We need not specify a concrete gauge group or its concrete representation though; the only assumption that we make is that the numbers of gluon and quark degrees of freedom are balanced in such a way that at zero temperature and density, the theory is in a confining phase where (approximate) chiral symmetry is spontaneously broken by the formation of a chiral condensate just like in QCD. As a consequence, the low-energy physics of the theory is dominated by the NG bosons of the global flavor symmetry.
It is well-known that for N flavors of massless quarks, (pseudo)real QCD-like theories possess an enlarged global symmetry, G = SU(2N ). This includes the usual chiral symmetry of QCD, SU(N ) L × SU(N ) R × U(1) B , as a subgroup. The additional symmetry transformations, only present in (pseudo)real theories, do not commute with the baryon number generator B, that is, convert quarks into antiquarks and vice versa. In the ground state, the symmetry is broken to a subgroup H ⊂ G by the formation of a chiral condensate. In real theories, H = SO(2N ), leading to 2N 2 + N − 1 NG bosons. In pseudoreal theories, on the other hand, H = Sp(2N ), leading to 2N 2 − N − 1 NG bosons. Out of these, N 2 − 1 are the usual pseudoscalar mesons, already present in QCD. The remaining NG modes correspond to scalar diquarks and antidiquarks, altogether N (N + 1)/2 pairs in real theories and N (N − 1)/2 ones in pseudoreal theories [23].
For the sake of an easy comparison with the phase structure of QCD in strong magnetic fields, we will restrict to (pseudo)real theories with N = 2 degenerate light quark flavors. (No further restrictions on the color degrees of freedom beyond those introduced above will be imposed.) It was shown in ref. [21] that such theories are free of the sign problem in the simultaneous presence of a baryon number chemical potential µ B and an external magnetic field H, provided that the electric charges of the u-type and d-type quarks satisfy the condition q u = −q d . Namely, the Euclidean Dirac operators D u,d of the two quark flavors are then related by the antiunitary mapping where C is the charge conjugation matrix, K the operator of complex conjugation, and P a matrix that realizes the similarity transformation between the generators of the gauge group and their complex conjugates. Eq. (2.1) implies that det D d = (det D u ) * , and thus ensures that the determinant of the full Dirac operator of the two-flavor theory is real and non-negative. The condition q u = −q d will be implicitly assumed in the following, unless explicitly stated otherwise. In general, the low-energy EFT for the NG degrees of freedom is constructed as a nonlinear sigma model on the coset space G/H. It would therefore appear from our above discussion that we need two different EFTs, one for real and one for pseudoreal theories. However, our task is dramatically simplified by recalling that we are interested in the possible presence of a CSL-like phase in the phase diagram, which requires a strong background magnetic field [1]. In a sufficiently strong magnetic field, electrically charged degrees of freedom become heavy due to Landau level quantization. The low-energy physics will then be dominated by the electrically neutral NG modes. To construct their EFT description, it is sufficient to consider the subgroups of G and H, left intact by the magnetic field. For two quark flavors and q u = −q d = 0, these are the same for real and pseudoreal theories, where U(1) Q refers to the subgroup, generated by the operator of electric charge Q. The low-energy spectrum contains accordingly dim G Q /H Q = 3 states, corresponding to the neutral pion and an electrically neutral quark-antidiquark pair. Note that the SU(2)×SU(2) subgroup of G Q is not the usual chiral symmetry of two-flavor QCD: it contains, among others, the operator of baryon number B as one of its "diagonal" generators.
To complete the setup for the construction of the low-energy EFT, we need to know the appropriate spacetime symmetry. In strong magnetic fields, assumed here, the Lorentz group SO(3,1) is broken explicitly down to its SO(1,1) × SO(2) subgroup; the former factor includes boosts along the magnetic field, while the latter rotations in the transverse plane. This reduced spacetime symmetry can be implemented by constructing a Lagrangian density in the usual manner, but instead of the Minkowski metric g µν contracting Lorentz indices using its projections, g µν to the two-dimensional subspace spanning time and the direction along the magnetic field, and g ⊥µν to the two-dimensional transverse plane.
Finally, we need a power-counting scheme to organize individual contributions to the effective Lagrangian. As usual in the chiral perturbation theory of QCD [24][25][26], we assign each derivative appearing in the Lagrangian density the order O(p 1 ). Likewise, the baryon chemical potential µ B counts as O(p 1 ). In contrast to the chiral perturbation theory, however, we will count the external magnetic field as O(p 0 ). This is necessary in order to keep consistency of the EFT in strong magnetic fields. Since the magnetic field is a singlet of both the internal symmetry G Q and the reduced spacetime symmetry, it is implemented by making the couplings in the effective Lagrangian arbitrary functions of H.
Having put together all the necessary pieces, we can now write down the effective Lagrangian at the leading, second order of the derivative expansion [27,28], Here Σ is the unitary unimodular 2 × 2 matrix field that includes the three degrees of freedom of the G Q /H Q coset space. Its covariant derivative defines the coupling of diquarks to the baryon number chemical potential, where b is the baryon number of a single quark and τ 3 the third Pauli matrix. The effective couplings in the Lagrangian (2.3)-the pion decay constant f π , pion mass m π and the velocity parameter v-are given by a priori unknown functions of the external magnetic field. Finally, the L WZ term represents the contribution of the chiral anomaly, known as the Wess-Zumino (WZ) term [29,30]. This has to be specified concretely before we can proceed with the analysis of the EFT.

Wess-Zumino term
In principle we could avoid having to construct the WZ term from scratch by using the results of ref. [31] and merely swapping the roles of baryon number B and electric charge Q therein. It is, however, instructive to outline at least the main steps. From the perspective of the low-energy EFT, we can at first forget about the anomaly in the underlying microscopic theory and simply ask the question whether there are nontrivial contributions to the effective Lagrangian that, unlike the first term in eq. (2.3), are invariant under the internal symmetry group G Q only up to a surface term [32,33]. It turns out that in four spacetime dimensions and for the internal symmetry (2.2), there is a single nontrivial WZ term, determined uniquely up to the addition of an arbitrary strictly invariant operator. It can be most conveniently expressed as a differential 4-form, ω WZ = A Q ∧ ω GW , where A Q is the gauge field 1-form of the U(1) Q factor of the symmetry group and ω GW is an invariant 3-form, corresponding to the so-called Goldstone-Wilczek current [34]. This is given explicitly for instance in eq. (4.5) of ref. [35], where F L,R are the field strength 2-forms associated with the left and right SU(2) factors of G Q , and the covariant derivative DΣ encodes the corresponding gauge field 1-forms, For our purposes, it is sufficient to consider gauge fields associated with the electric charge Q and baryon number B, A straightforward manipulation then leads to an expression for the WZ Lagrangian, where C is an as yet undetermined normalization factor.
The latter can be fixed by considering a single operator contributing to the WZ term (2.7) and matching it to the underlying gauge theory. A minor modification of the standard derivation of the Abelian anomaly (see section 22.2 of ref. [36]) gives a contribution to the effective Lagrangian, linear in the neutral pion field π 0 , as where d is the dimension of the representation of the color gauge group that a single quark transforms in, and the remaining trace is performed over the space of quark flavors. The field strength F µν is diagonal in flavor and includes the contributions of both the electric charge and the baryon number. For future reference, we also include a contribution from an external gauge field coupled to isospin, that is, write schematically F µν = QF Q µν + BF B µν +IF I µν . Inserting finally the proper values of the quantum numbers, This result applies to all QCD-like theories with two quark flavors regardless of the color gauge group; setting d = 3, b = 1/3, q u = 2/3 and q d = −1/3 reproduces the known result for QCD itself [31]. To match it to the form of WZ term (2.7) valid for (pseudo)real QCDlike theories, we insert Σ = exp i fπ τ 3 π 0 therein and integrate by parts. For q u = −q d as assumed in this paper, the purely electromagnetic operator π 0 µναβ F Q µν F Q αβ , responsible for the two-photon decay of neutral pions, is missing. Comparing the coefficients of the operator µναβ F Q µν A B α ∂ β π 0 in eqs. (2.7) and (2.9) then gives The construction of our EFT is completed by replacing A B µ with (µ B , 0). Note that due to the modified power counting for the background magnetic field, the WZ term (2.7) indeed belongs to the leading order Lagrangian (2.3) as advertised.

Effective Lagrangian and Hamiltonian
The remainder of the paper is devoted to an analysis of the EFT as defined by eqs. (2.3), (2.7) and (2.10). We will from now on assume without loss of generality that q u = −q d = 1/2 and b = 1/2. This convention is natural in that it maintains the physical electric charges of the three pions, and gives the diquarks a unit baryon number. Any other choice can be absorbed into a redefinition of H and µ B , respectively.
To prepare for the analysis, we first rewrite the Lagrangian (2.3) in terms of a unit vector variable (n 0 , n), defined by This brings the Lagrangian to the form where we oriented the magnetic field along the z-axis, and subtracted a constant to ensure that the Lagrangian vanishes in the vacuum where (n 0 , n) = (1, 0). Finally, the ellipsis stands for three-derivative contributions coming from the WZ term (2.7). We will drop these contributions, as they have no effect on either the structure of the ground state or the excitation spectrum. Indeed, choosing a gauge where A Q 0 = 0, it is easy to see that they are necessarily linear in time derivatives, and thus do not contribute to the Hamiltonian. Moreover, we will show below that the ground state of the EFT is either uniform or modulated in the direction of the magnetic field. When analyzing the excitation spectrum, each of the three derivatives in the three-derivative terms then either acts on the z-dependence of the ground state, or generates a factor of momentum; their contributions then vanish by antisymmetry.
For the record, we also write down another form of the Lagrangian, where the unit vector (n 0 , n) is mapped on three angular variables via n 0 = cos θ cos φ, n 1 = sin θ cos α, n 2 = sin θ sin α, This converts the Lagrangian (2.12), now dropping the three-derivative terms, into (2.14) Upon carrying out a Legendre transform, it is easy to see that derivatives with respect to time and the transverse directions only enter the Hamiltonian in squares. The ground state is thus necessarily independent of both time and the transverse coordinates, and can be found by minimization of the spatial average of the one-dimensional effective Hamiltonian, (2.15) 3 Chiral limit

Ground state
In the chiral limit where m π = 0, the ground state can be found straightforwardly by rewriting the Hamiltonian (2.15) as where we introduced the shorthand notation From the first line of eq. (3.1), it is obvious that the ground state is realized for constant θ. At the same time, the expectation value of the angle φ should satisfy ∂ z φ =Hµ B as long as cos θ = 0. The second line of eq. (3.1) picks the most favorable value of θ, depending on whetherH is smaller or greater than one. Altogether, we find two different phases: •H < 1: the ground state is at θ = π/2, corresponding to n 0 = 0, n 1 = cos α , n 2 = sin α , with constant α . This describes a Bose-Einstein condensate (BEC) of diquarks.
•H > 1: the ground state is at θ = 0, corresponding to up to a shift of the z-coordinate. This describes the CSL state in the chiral limit [1].
At any nonzero value of the baryon chemical potential, the phase diagram of the EFT in thē H-µ B plane (and at zero temperature) therefore contains two phases. The µ B -independent critical magnetic field for the formation of the inhomogeneous CSL state is, by eq. (3.2), The parametric dependence on the number of quark color degrees of freedom, encoded in the constant d, arises from the competition of the anomalous and non-anomalous parts of the effective Hamiltonian, only the former of which is proportional to d. We will return to this observation in the conclusions.

Excitation spectrum
The angular parameterization (2.13) becomes singular at both θ = 0 and θ = π/2. In order to discuss the fluctuations about the ground state, we therefore have to return to the unit-vector parameterization, in which the Lagrangian takes the form (2.12). Let us first focus on the BEC phase (3.3). By setting α without loss of generality to zero, we can parameterize n 1 and n 2 as The effective Lagrangian (2.12) can then be expanded in powers of the independent fluctuations n 0 , n 3 and α. Keeping only terms bilinear in these fluctuations, we obtain Obviously, α describes a gapless pseudo-relativistic excitation with the dispersion relation This is the NG boson of the spontaneously broken U(1) B symmetry. The n 0 and n 3 modes are mixed by the anomaly term. Diagonalizing the 2 × 2 inverse propagator leads to a pair of dispersion relations, Both these modes become gapless at the phase transition to the CSL phase, whereH = 1.
The nonuniform nature of the CSL state is manifested by the fact that the minimum of ω ± (p) appears at a nonzero value of p z , exactly corresponding to the momentum scale of the CSL solution (3.4). Note that the energy of these excitations at zero momentum is fully determined by the baryon chemical potential alone, ω ± (p = 0) = µ B , in accord with the interpretation of both modes as massive NG bosons [37].
Moving on to the CSL phase, we need to find a suitable parameterization of fluctuations about the ground state (3.4). This time, we can set The three independent excitation modes then correspond to n 1 , n 2 and φ. The bilinear part of the Lagrangian (2.12) is then, upon some manipulation, seen to be We find again one gapless pseudo-relativistic excitation with the dispersion relation This pion-like mode can be interpreted as the NG boson of the spontaneously broken translation invariance in the CSL state. The n 1 and n 2 modes describe the diquarkantidiquark pair, and their dispersion relation is obtained by diagonalization, We can see that the diquark mode becomes gapless at the phase transition,H = 1, signaling an instability towards the BEC state.

Away from the chiral limit
In order to investigate the phase diagram away from the chiral limit, we take advantage of the angular parameterization (2.13). The analysis is simplified by introducing dimensionless spacetime coordinatesx µ and a dimensionless chemical potential x via 2 The effective one-dimensional Hamiltonian (2.15) then also acquires a dimensionless form, where the primes denote derivatives with respect toz.

Ground state
Finding the absolute minimum of (the spatial average of)H seems difficult in a closed form, we will therefore have to make some further assumptions about the structure of the ground state. We have performed extensive checks, both numeric and where possible analytic, that no state with energy lower than the states discussed in detail below exists for any choice of the dimensionless parameters x andH. Some details of these checks are provided for an interested reader in appendix A. The analysis becomes elementary if we for a moment restrict to uniform field configurations. Finding the ground state then boils down to a minimization of the static part of the Hamiltonian (4.2) with respect to θ and φ. We find two different stationary states: The BEC state only exists as a stationary state of the Hamiltonian if x > 1, or µ B > m π , and it automatically has a lower energy than the trivial vacuum. This exhausts the possible candidate ground states in the subspace of uniform field configurations. To prove that, for given values of x andH, the ground state is inhomogeneous, it is therefore sufficient to find any state that has a lower energy than the two states listed above.
To that end, we inspect a different class of states by setting θ → 0. Upon this replacement, the effective Hamiltonian (4.2) reduces to that of QCD where only the neutral pion degree of freedom survives in the low-energy EFT. The minimum of the Hamiltonian on such a subset of field configurations is thus known, and corresponds to the CSL state [1], where sn is one of Jacobi's elliptic functions and k the associated elliptic modulus. This is fixed by a minimization of the spatially averaged Hamiltonian, which leads to the condition, where E(k) is the complete elliptic integral of the second kind. Using the explicit solution (4.5) along with eq. (4.6), the spatially averaged energy density of the CSL state can be cast asH Given that the elliptic modulus falls into the range 0 ≤ k ≤ 1, this state always has a lower energy than the trivial vacuum. Since E(k)/k ≥ 1, it only exists forH ≥ 4/(πx) though.
Comparing the energies of all the three candidate states gives rise to the phase structure shown in figure 1. We can see that there is a range of magnetic fields where, as anticipated, the ground state becomes inhomogeneous. The critical magnetic field for the formation of an inhomogeneous ground state is given by (4.8) Note that the limit x → ∞ is equivalent to the chiral limit, m π → 0. For large x, the critical magnetic field separating the BEC and CSL phases therefore goes toH = 1.
Except for the vacuum phase, the ground state is characterized by nonzero baryon number density, which can be obtained from the Hamiltonian (2.15) by taking a derivative with respect to the chemical potential. Defining a rescaled, dimensionless baryon density asn B ≡ n B /(f 2 π m π ), its value in the BEC phase, and the spatially averaged value in the CSL phase [1], can be expressed as 3 where K(k) is the complete elliptic integral of the first kind, and k is determined implicitly by eq. (4.6). The dependence of the baryon number density on magnetic field for several different values of chemical potential is shown in figure 2. The numerical results show that while in the BEC phase, the baryon density is driven solely by the chemical potential, in the CSL phase, it is further boosted by increasing magnetic field. This is in accord with eq. (3.1) which implies that in the chiral limit, n B = f 2 π µ BH 2 . The corresponding dimensionless density,n B = xH 2 , describes accurately the curves shown in figure 2.

Excitation spectrum
It is of course possible that by restricting to the three candidate states as discussed above, we miss a phase in the phase diagram where a more sophisticated, necessarily inhomogeneous, ground state has an even lower energy. Apart from the direct checks reported in appendix A, one can also check independently the consistency of the phase diagram displayed in figure 1 by evaluating the spectrum of excitations above the tentative ground state. Should this reveal an instability in some of the dispersion branches, we would have a solid evidence for the existence of a new, more favorable ground state. Let us start with the region in the phase diagram where the trivial vacuum seems to prevail. This corresponds to n 0 = 1 and n = 0. To extract the spectrum, we return to eq. (2.12) and use n as the three independent degrees of freedom along with n 0 = √ 1 − n 2 .

The bilinear part of the effective Lagrangian (2.12) then reads
L vac bilin = f 2 π 1 2 (g µν +v 2 g µν ⊥ )∂ µ n·∂ ν n+µ B (n 1 ∂ 0 n 2 −n 2 ∂ 0 n 1 )+ This is a known result: the spectrum contains a neutral pion with mass m π and a diquarkantidiquark pair with mass m π and chemical potential µ B . We can readily identify the instability towards diquark BEC at µ B = m π . On the contrary, there is nothing in the spectrum that would indicate the onset of the inhomogeneous CSL phase.

Spectrum in the BEC phase
The spectrum of the BEC phase can be analyzed in the same manner, except that we can now take advantage of the angular parameterization (2.13). Denoting the value of θ in the BEC ground state as θ 0 and expanding in the fluctuations of all the three angles θ, φ, α, a straightforward manipulation based on eq. (2.14) gives the rescaled bilinear Lagrangian

(4.11)
This readily leads to the inverse matrix propagator in the (θ, φ, α) space, where we used the shorthand notation¯ v ≡ (g µν +v 2 g µν ⊥ )∂ µ∂ν . Upon Fourier transforming to the space of frequencyω and momentump, this becomes where we likewise abbreviatedp 2 v ≡ v 2p2 ⊥ +p 2 z . The dispersion relations of the three excitation branches in dimensionless form can in principle be obtained by setting the determinant of the inverse propagator in eq. (4.13) to zero. It is unfortunately not possible to write the result, corresponding to the solution of a cubic equation forω 2 , in a closed real form. 4 However, what interests us most is whether some of the modes becomes gapless in a part of the tentative BEC phase. To search for such an instability, it is sufficient to setω = 0, upon which the determinant of the inverse propagator is easily evaluated and factorized, 14) The leading factor ofp 2 v reflects the presence of the gapless mode due to the spontaneously broken exact U(1) B symmetry. The expression in the square brackets is manifestly positive forH = 0, we are thus looking for the lowest value ofH,H inst , such that drops to zero for somep. A straightforward function analysis leads to the result, The value ofH inst is indicated by the upper dashed line in figure 1. This confirms that the diquark BEC state is stable under small field fluctuations everywhere in the tentative BEC phase. The critical fieldH inst merely denotes a point where the BEC state ceases to be a local minimum of the Hamiltonian inside the CSL phase; this defines the spinodal curve for the first-order phase transition between the BEC and CSL phases. In order to illustrate the onset of instability as the magnetic field approachesH inst , we show in figure 3 the evolution of the dispersion relation of the lightest mode in the spectrum, obtained by numerically solving the condition that the determinant of the inverse propagator (4.13) is zero. We can see that close to the critical magnetic field, the dispersion relation develops a characteristic roton-like minimum, which acts as a precursor to the instability with respect to condensation of modes with nonzero momentum.
A cautious reader might be puzzled by the fact that the instability atH inst occurs in the dispersion relation of the lightest mode, which corresponds to the NG boson of the  spontaneously broken exact U(1) B symmetry. Namely, our previous analysis of the chiral limit in section 3.2 clearly showed that the instability appears in the dispersion relation of modes that have a nonzero gap at p = 0, not of the NG mode, see eqs. (3.8) and (3.9). The resolution of this apparent contradiction lies in the mixing of all three modes, present away from the chiral limit. This is illustrated by the upper two panels of figure 4, where we display the dispersion relations of all three light modes for fixed x and two different values ofH. The results in the chiral limit are indicated by the dashed lines. 5 We can see that switching on a nonzero quark mass leads to avoided level crossing in the spectrum, as a result of which the roton-like minimum indeed appears in the dispersion relation of the lightest mode.
To get further insight into the nature of the mixing, we extracted the coupling of the lightest mode to the θ, φ and α fields using the Källén-Lehmann spectral representation. Schematically, the inverse of eq. (4.13) near a poleω = (p) takes the form where χ i runs over θ, φ, α. In the lower two panels of figure 4, we show the values of c i ∝ | 0|χ i |p | 2 , normalized so that c θ + c φ + c α = 1. With this definition, the displayed coefficients highlight the composition of the lightest mode in the spectrum. Obviously, it is "mostly α" before the avoided crossing (in accordance with the NG nature of the mode), and "mostly φ" after the avoided crossing.

Spectrum in the CSL phase
The spectrum in the CSL phase can be analyzed analogously, except that we now have to deal with the fact that we are looking for fluctuations of an inhomogeneous ground state. Given that θ = 0 in the CSL state (4.5), the parameterization (2.13) is singular and we have to return once more to the unit vector (n 0 , n). Denoting the spatially-dependent expectation value of φ in the CSL ground state as φ 0 , we can use a parameterization that straightforwardly generalizes eq. (3.10), The three independent degrees of freedom can be taken as n 1 , n 2 and φ. However, since the exact baryon number U(1) B symmetry remains unbroken in the CSL state, it is more convenient to trade n 1 , n 2 for the baryon number eigenstates, The bilinear part of the effective Lagrangian then splits into separate pieces for φ and n ± , L CSL bilin = L CSL φ + L CSL n ± . Upon rescaling, the two pieces take the form The bilinear Lagrangian L CSL φ of the neutral pion sector is, apart from the appearance of the velocity parameter v, identical to that previously found in QCD [1], and we can therefore merely quote the result for the spectrum. It consists of two energy bands, the lower of which is gapless in accord with the fact that the CSL state spontaneously breaks translations in the direction of the magnetic field. In the gapless (valence) band, the transverse motion is pseudo-relativistic with phase velocity v. The dispersion relation is, however, nonlinear in the longitudinal direction, that is, along the magnetic field. The phase velocity at long wavelengths is given by eq. (5.3) of ref. [1]. Obviously, no instability is present in the neutral pion sector.
The L CSL n ± part of the bilinear Lagrangian describes the propagation of (anti)diquarks. Upon Fourier transform in time and the transverse coordinates, one readily finds that their (dimensionless) dispersion relations can be expressed as where λ runs over the eigenvalues of the differential operator Upon a further rescaling of the coordinate,z ≡z/k, using eqs. (4.5) and (4.6) and the properties of the Jacobi elliptic functions, the eigenvalue problem for ∆ can be recast as whereλ runs over the eigenvalues of another operator, dn being another of the Jacobi elliptic functions.
In order to see whether there is a region in the CSL phase where the diquark mode becomes unstable, we have to check whether for some values of x andH, the smallest eigenvalue λ drops down to x 2 . The operator∆, and thus also the eigenvalues λ, depends only on the elliptic modulus k. The most convenient way to find the condition for the (in)stability of the diquark on the CSL background is therefore to fix the value of x, find the k for which the smallest eigenvalue of ∆ equals x 2 , and then determine the magnetic field from eq. (4.6).
It remains to find the ground state, that is the lowest eigenvalue, of∆ as a function of k. This can be done numerically using the variational principle, and the results are shown in figure 5. Some analytic insight into the results can be achieved by looking at the asymptotic behavior of the potential V (z, k) as a function of k. First, we have lim k→0+ V (z, k) = 8. Thus, at k = 0, the operator∆ describes a free particle with spectrum starting at minλ = 8. Second, we find that Remarkably, this potential does not have any bound state in spite of its minimum at the origin. 6 Accordingly, it has a single energy band, starting atλ = 6. We thus conclude that lim k→0+ (minλ) = 8 and lim k→1− (minλ) = 6. 6 This was verified by a numerical examination of the spectrum of the one-parameter set of potentials Vα(z) ≡ 6 tanh 2z + 16α/(π coshz). For α = 0, this is equivalent to the well-known Pöschl-Teller potential, whose spectrum can be found analytically. It consists of two bound states with eigenvalues 2 and 5, and a continuous spectrum above 6. For 0 < α < 1, the spectrum was investigated numerically using the variational method. With increasing α, the eigenvalues corresponding to the bound states grow and eventually disappear in the continuum: the upper bound state disappears for α ≈ 0.4, the lower bound state later at α ≈ 0.9. Hence for α = 1, corresponding to the potential (4.25), there is no bound state left. The analysis for small k can be further refined by considering the first nontrivial order in the series expansion of the potential in powers of k, Using the numerical results displayed in figure 5 leads to the stability limit for diquarks on the CSL background, indicated by the lower dashed line in figure 1. We observe that the diquark mode only becomes unstable at magnetic fields below the transition to the BEC phase. In other words, the CSL phase as indicated in figure 1 is stable under small fluctuations; the lower dashed line merely denotes the spinodal curve where the CSL state ceases to be a local minimum of the Hamiltonian.
The dispersion relation (4.21) combined with the data shown in figure 5 can be used without delving further into the spectrum of the operator ∆ to determine the mass spectrum of diquarks. Indeed, by eq. (4.21), the (rescaled) mass, that is energy at zero momentum, of the diquark-antiquark pair equals min λ ∓ x. In figure 6, we show the mass of the lighter of the two states as a function of magnetic field for several fixed values of the baryon chemical potential. This further highlights the location of the transition between the BEC and CSL phases. For higher x, the discontinuity is so small that it is not visible in the graph.

Summary and conclusions
In this paper, we have analyzed in detail the effect of chiral anomaly on the phase diagram of (pseudo)real QCD-like theories subject to a strong external magnetic field. We showed that under fairly general assumptions, their phase diagram features a nonuniform phase of the CSL type, where nonzero baryon density is generated by a topological crystalline condensate of neutral pions, similarly to what was previously observed in QCD [1]. All the important results were discussed in detail in the main text of the paper, here we therefore merely append a few concluding remarks. First, as already explained in the prequel to this paper [21], the main result as displayed in figure 1 should be taken with a grain of salt. Namely, all the parameters of the effective Lagrangian (2.3) are given by a priori unknown functions of the magnetic field. Given that, for instance,H = CH/[f π (H)] 2 by eq. (3.2), the variables used to label the axes in figure 1 represent a complicated nonlinear mapping of the physical parameters µ B and H. It is therefore not a priori obvious whether the regionH > 1, where the inhomogeneous CSL phase appears, can in fact be reached for any finite value of H.
Based on eq. (3.5), one can argue that in theories with sufficiently large gauge group and its sufficiently large representation (that is large d), the critical magnetic field for the transition to the CSL phase is small and thus under theoretical control. Indeed, the conditionH > 1 can be rewritten as H/(4πf π ) 2 > 1/(2d). Combined with the fact that 4πf π is the loop factor that controls the derivative expansion of the EFT [39], this suggests that the critical field for the formation of the CSL is within the reach of the EFT. Using an explicit one-loop expression for the H-dependence of f π , available in the literature, shows that f π can, in fact, be treated as a constant; the correction to the critical magnetic field ensuing from its H-dependence is suppressed by 1/d and thus negligible for large d [21].
For theories with small d, most notably for two-color QCD, our results are somewhat inconclusive though. In order to decide whether theH > 1 region can be reached with a finite magnetic field, one would need detailed input on the magnetic field dependence of the pion decay constant f π (H). Such an input is, to the best of our knowledge, currently not available.
Second, by focusing solely on theories free of the sign problem, we made a prediction for the phase diagram of a class of QCD-like theories that can serve as a reference for future lattice simulations of inhomogeneous phases of dense quark matter. On the lattice, the presence of the CSL phase can in principle be tested in several different, more or less explicit ways. The most straightforward, but also most brute force, approach would be to look for the nonuniform order, either by trying to study directly the order parameter, or its correlator at long distances [19]. Alternatively, one could supplement the EFT analysis presented here by a nonperturbative evaluation of the neutral pion decay constant as explained above. Last but not least, the presence of the characteristic roton-like mode in the spectrum, as displayed in figures 3 and 4, makes it possible to demonstrate the instability towards the formation of an inhomogeneous ground state by studying detailed properties of the adjacent homogeneous BEC phase. 7 Finally, although we did focus solely on theories free of the sign problem, CSL is clearly a much more general phenomenon, only requiring an anomalous coupling of an electrically neutral degree of freedom to an external magnetic field and a chemical potential. Apart from QCD at nonzero baryon chemical potential studied in ref. [1], this applies in particular to QCD at nonzero isospin chemical potential µ I . This theory is known to be free of the sign problem in a vanishing magnetic field [40], 8 and its phase diagram has been studied thoroughly with an increasing precision, using both lattice simulations [42] and EFT [43]. By adding a gauge field for isospin to the derivation of the WZ term in section 2.1, one can see that an external magnetic field in combination with the isospin chemical potential can also give rise to the CSL phase. The resulting low-energy EFT is identical to that analyzed in ref. [1] except for a different overall normalization of the WZ term, fixed by eq. (2.9).
A similar comment applies to the class of QCD-like theories covered here but with nonzero quark electric charges that do not satisfy the condition q u = −q d , assumed throughout our paper. In this case, all the diquarks are electrically charged and thus acquire a gap from the magnetic field. The low-energy spectrum then only contains the neutral pion, and the analysis of ref. [1] applies without modification, possibly only with a different coefficient of the WZ term. 7 The utility of the roton-like dispersion relation as a smoking gun that hints at a presence of a nonuniform phase was pointed out to us by Gergely Endrődi. 8 Switching on a magnetic field destroys the positivity of the determinant of the Dirac operator unless one tunes the quark electric charges so that qu = q d [41]. In this case, all three pions are electrically neutral and the resulting low-energy EFT is equivalent to that studied in this paper up to a replacement of the (anti)diquark degrees of freedom with the "charged" pions π ± and of µB with µI . The only difference is in the value of the coefficient C; here eq. (2.10) should be replaced with C = d 8π 2 (qu + q d ).
where the function k(θ) is defined implicitly by a condition, generalizing eq. (4.6), This equation has a unique solution for k only for cos θ > (4/πxH) 2 . For values of θ not satisfying this condition, k = 1 is to be used in the Hamiltonian (A.4). It remains to minimizeH(θ) as defined by eq. (A.4) with respect to θ for given values of x andH. As it turns out, however, the ground state has, for any x andH, either φ = 0 (BEC state) or θ = 0 (CSL state). Within the class of states considered here, there is no other state of even lower energy.
A.2 States with periodic θ and linear φ Our first ansatz discussed above is probably not the most natural one. After all, once translational invariance is spontaneously broken by a spatially varying φ as in the CSL state, one would expect θ, if nonzero, to be nonuniform as well. In this case, we have to make some further assumptions on the spatial profile of both θ and φ in order to be able to compute, and minimize, the energy quasi-analytically.
To that end, we first recall the physical meaning of the angular variables θ and φ as introduced by eq. (2.13). Since we can set α = 0 without loss of generality, θ and φ are nothing but the usual spherical coordinates parameterizing a unit two-sphere in the n 0 -n 1 -n 3 space. The angle θ measures the azimuthal deviation from the equator, lying in the n 0 -n 3 plane, whereas φ is the usual polar angle. The BEC state maps to a fixed point in the n 0 -n 1 plane, whereas the CSL state describes a periodic motion along the equator with the z-coordinate playing the role of "time." A nonuniform state that should have a lower energy than the uniform BEC state, must make use of the anomalous contribution to the Hamiltonian (4.2), proportional to φ , and thus have nontrivial topology due to winding in the φ-direction. We want to see if it is possible to lower the energy as compared to the CSL state by allowing the azimuthal angle θ to vary as well. Due to the nonperiodic nature of the azimuthal angle, we expect θ to oscillate within some finite range.
The simplest ansatz with these properties takes the form where θ 0 , φ 0 , K and l are parameters. The linear ansatz for φ is motivated by the profile of the CSL state in the chiral limit (3.4), which is combined with the simplest ansatz for θ oscillating with the same period as the CSL state. Note that this ansatz contains neither the BEC nor the CSL state as its special cases. First, the spatial profile of φ(z) is fixed to a simple approximation to the more complicated CSL state. Second, the ansatz (A.6) has a nontrivial topology for any finite l; the BEC state can be formally recovered by setting l = ∞, which however does not correspond to the limit of eq. (A.6) as l → ∞.