Influence of vector interactions on the favored shape of inhomogeneous chiral condensates

We update a previous study of the effects of vector interactions in the Nambu--Jona-Lasinio model on the formation of inhomogeneous chiral symmetry breaking condensates. In particular, by properly considering a spatially modulated vector mean-field associated with the quark number density of the system we show that, as the value of the vector coupling increases, a chiral density wave modulation can become thermodynamically favored over a real sinusoidal modulation. This behavior is found both via a Ginzburg-Landau analysis close to the Lifshitz point, as well as with a full numerical diagonalization of the mean-field Dirac Hamiltonian at vanishing temperature.


I. INTRODUCTION
The phase structure of QCD at high densities and intermediate temperatures is expected to be extremely rich. A particularly fascinating scenario suggested by explicit model calculations revolves around the formation of an inhomogeneous "island", characterized by the presence of crystalline chiral condensates, which can appear in the region where the firstorder chiral phase transition is expected to take place (for a recent review, see [1]). As a consequence, the critical point (CP) as a cornerstone of the phase diagram is replaced by a Lifshitz point (LP) where three different phases -a homogeneous broken, an inhomogeneous and a chirally restored one -coexist. Among various shapes which have been compared so far, the energetically favored crystalline structure within the inhomogeneous island is the so-called "real-kink crystal" (RKC) [2,3], which assumes the form of an array of domain-wall solitons close to the phase transition to the homogeneous phase and rapidly evolves into a sinusoidal shape as one approaches chiral restoration. The existence of an inhomogeneous phase might have phenomenological consequences for compact stars [4] and for the physics of heavy-ion collisions, especially for the future facilities at FAIR and NICA as well as the beam-energy scan at RHIC [5,6].
A common extension considered in effective models of QCD such as the Nambu-Jona-Lasinio (NJL) model for describing matter at high densities is the inclusion of a repulsive vector interaction channel [7][8][9][10][11]. Within the mean-field approximation, the time component of the vector condensate can be identified with the density of the system, and induces an effective shift in the chemical potential of the quarks [7]. If only homogeneous matter is considered, aside from moving the chiral transition line to higher chemical potentials, the most notable effect of vector interactions is a shift of the critical point towards lower temperatures, leading to its disappearance beyond a certain value of the vector coupling constant [10].
A first investigation of the effects of vector interactions on chiral crystalline phases has been performed in Ref. [12]. There it was found that the inclusion of vector interactions leads to a dramatic enhancement of the inhomogeneous region in the model phase diagram and the disappearance of the critical point inside the inhomogeneous phase. In that study, in order to be able to use known analytical expressions for the eigenvalue spectrum of the Dirac Hamiltonian in the presence of certain crystalline backgrounds, instead of considering a spatially modulated density, the vector mean field was approximated by its spatial average. While this approach is still self-consistent when a plane-wave modulation of the chiral condensate (typically referred to as a "(dual) chiral density wave" (CDW) [13]) is considered, this is certainly not the case for modulations like the RKC [12] or even a simple sinusoidal ansatz. For these cases, a repulsive vector interaction might generate an additional energy cost for the formation of a spatially modulated quark density, possibly influencing the competition between different crystalline phases and altering the resulting phase structure of the model. In this work, we aim at investigating this possibility by lifting the assumption of constant density and allow for a more consistent treatment of spatially modulated vector condensates within the model. This paper is structured as follows: We begin by describing our model and the procedure for including inhomogeneous chiral and vector condensates in Sec. II, while in Sec. III and Sec. IV we present results from a Ginzburg-Landau investigation and from a numerical diagonalization of the Dirac Hamiltonian, respectively. We finally conclude in Sec. V.

II. MODEL AND BASIC FORMALISM
We consider an NJL model in the chiral limit, given by the Lagrangian where ψ denotes the quark field with vanishing bare mass, N f = 2 flavor and N c = 3 color degrees of freedom. In addition to the standard isoscalar scalar and isovector pseudoscalar interactions with coupling constant G S , the Lagrangian contains an interaction in the isoscalar vector channel with coupling constant G V . Here γ µ are Dirac matrices, and τ a are Pauli matrices in flavor space.
Performing the mean-field approximation, we allow for the scalar ψ (x)ψ(x) ≡ S(x), pseudoscalar ψ (x)iγ 5 τ a ψ(x) ≡ P (x)δ a3 and vector ψ (x)γ µ ψ(x) ≡ n(x)g µ0 condensates, which are restricted to be time-independent but allowed to vary in space in order to describe possible inhomogeneities. As in Ref. [12], we assume that the pseudoscalar condensate is diagonal in flavor space and consider only the 0-component of the vector condensate, corresponding to the quark-number density n(x) = ψ (x)γ 0 ψ(x) . This latter assumption is well justified when considering spatially uniform matter since non-vanishing spatial vector components break the rotational invariance of the system. While this argument is no longer valid in phases where the rotational symmetry is already broken by the scalar condensates, vanishing space-like components of the vector condensate could also be protected if the ground state is invariant under parity transformations. This is the case for all modulations considered in this article. Indeed preliminary studies seem to confirm that ψ γ k ψ = 0, k = 1, 2, 3, holds for these phases [14,15].
Following standard procedures (see Ref. [1] for details), the mean-field thermodynamic potential per unit volume V at a given temperature T and chemical potential µ can then be written as a sum of a kinetic and a condensate contribution, where the latter is simply given by a volume integral over the condensates, The kinetic contribution, basically amounts to a sum over the eigenvalues E λ of the mean-field Dirac Hamiltonian where corresponds to the mean-field Hamiltonian without vector interactions. Here is a complex constituent quark mass function and P ± = 1 2 (1 ± γ 5 τ 3 ) are chirality projectors. Generalizing a method often employed in homogeneous phases [7], we can also absorb the effects of the vector mean field in Eq. (5) into a space-dependent shifted chemical potential function,μ so that The condensate part of the thermodynamic potential, Eq. (3), is then given by Already for the simpler case without vector interactions, the diagonalization of the Dirac Hamiltonian in the presence of inhomogeneous condensates is a non-trivial task, as an explicit space dependence of the mean fields leads to a non-diagonal structure in momentum space.
For two particular one-dimensional spatial modulations of the order parameters, however, it is possible to obtain an analytical expression for the eigenvalue spectrum of the Hamiltonian when vector interactions are switched off. The first one is the CDW [13], which is a simple plane-wave modulation where, without loss of generality, we have chosen the wave to vary along the z-direction.
The amplitude ∆/ √ 2 and the wave number q have to be detemined as functions of T and µ by minimizing the thermodynamic potential.
The other example is the one-dimensional RKC, which can be expressed in terms of Jacobi elliptic functions, with parameters ∆ and ν, which are again determined by minimizing Ω. The so-called elliptic modulus ν thereby determines the shape of the mass function. It is found that this type of solution assumes a solitonic shape close to the onset of the inhomogeneous phase at low chemical potentials, while as the chemical potential increases it quickly assumes a sinusoidal shape, see [12,16].
In the absence of vector interactions, the RKC solution is energetically favored over the CDW throughout the entire inhomogeneous window [2]. It is also more favored than all other shapes (including higher-dimensional modulations) which have been tested so far either by numerical diagonalization of the Hamiltonian [3] or within Ginzburg-Landau studies [17,18].
In general, the inclusion of vector interactions alters the energy spectrum of the model through the appearance of the new mean-field n(x) in the Hamiltonian, and therefore the known analytical expressions for the eigenvalues of H 0 cannot be used. A notable exception however occurs if the density of the system is spatially homogeneous: In this case, as obvious from Eq. (5), the eigenvalues of H 0 are just shifted by a constant amount 2G V n, and the only additional task is to determine the value of n selfconsistently. In fact, for constant density the shifted chemical potentialμ is also constant, so that Ω kin (T, µ) in the presence of vector interactions is identical to Ω kin (T,μ), evaluated for G V = 0. The selfconsistent value ofμ (and thus of n) is then obtained from the stationarity of the thermodynamic potential, ∂Ω/∂μ = 0.
However, while for the CDW ansatz n(x) is indeed constant, this is not the case for the RKC. In Ref. [12], in order to nevertheless make use of the known analytical expressions for the eigenvalue spectra, the spatially modulated vector condensate in Eq. (5) was approximated by its spatial averagen ≡ n(z) . Hence, within this "average-density approximation" (ADA), the above feature that the main effect of the inclusion of vector interactions is just a shift µ →μ in Ω kin holds for all modulations, including the RKC.
In turn, this means that within the ADA the competition between different types of modulations is left substantially unaltered compared to the G V = 0 case. In particular the RKC remains favored over the CDW [12]. On the other hand, since the true density in a RKC is not constant, while it is in a CDW, one may ask to what extent this result is just an artifact of the approximation. In fact, since the vector interaction is repulsive, it is conceivable that homogeneous density distributions, as in the CDW, eventually become favored over modulated ones, as in the RKC, when the vector coupling is increased. In the following we will investigate this possibility.

III. GINZBURG-LANDAU EXPANSION
In the vicinity of the LP, where the homogeneous and inhomogeneous phases with broken chiral symmetry meet with the restored phase, the mass function M (x) and its gradient are small, so that some insight into the properties of the system can be inferred via a Ginzburg-Landau (GL) analysis of the thermodynamic potential. In this section we want to apply this technique to study the effect of vector interactions.
To that end, following [12], we write the shifted chemical potential as the sum of its value in the restored phase and a small deviation,μ(x) =μ rest 0 + δμ(x), and expand Ω in powers of the mass function M (x) and δμ(x) around the chirally restored solutions (M,μ) = (0,μ rest 0 ). More specifically, we write Ω = Ω 0 + 1 V d 3 x Ω GL with Ω 0 the thermodynamic potential in the restored phase and where the terms with the coefficients α i (T, µ) correspond to the potential for G V = 0, while the terms with the coefficients β i (T, µ) appear when vector interactions are included. A remarkable feature of the former part is that the coefficients α 41 of |M | 4 and α 42 of |∇M | 2 are equal (≡ α 4 ), with the consequence that for G V = 0 the LP, given by α 2 = α 42 = 0, coincides with the CP, which is given by α 2 = α 41 = 0 [19].
In Eq. (13) we included terms up to order M 6 in order to be able to determine the favored type of modulation for the chiral condensate [17,19]. Our power counting is given by assuming that ∇ ∼ M , as customary, while δµ is of order M 2 [12]. Indeed, by solving the Euler-Lagrange equation for δμ we find to up to order M 4 where we defined R β = β 42 /β 41 . To leading order we thus have δμ ∼ M 2 , as claimed above.
Note that for a CDW modulation, Eq. (11), all the quantities on the right-hand side of Eq. (14) are spatially constant, so that the shifted chemical potential and therefore the density are constant as well.
By plugging the solution of Eq. (14) into Eq. (13) and truncating again at O(M 6 ) we arrive at The computation of the GL coefficients is tedious but straightforward. The α i have already been obtained in Ref. [19] for the case G V = 0. Taking into account the shift µ →μ rest 0 due to the vector condensates one finds for G V = 0 and thus which is negative for positive values of G V . Note that these coefficients, which multiply powers of δμ, do not have any divergent parts and therefore remain unregularized. This can easily be understood by observing that these terms can be obtained from the thermodynamic potential for homogeneous order parameters as derivatives with respect toμ. Hence they originate from the unregularized medium contribution to the thermodynamic potential and therefore remain unregularized as well. This argument does not apply to gradient terms, which can in principle carry a regularization dependence.
From Eq. (15) we can now see that when G V (and consequently R β ) is nonzero, the GL coefficient multiplying |M | 4 acquires an additional contribution proportional to R β which leads to the known shift of the CP to lower temperatures [10,12]. On the other hand, the |∇M | 2 coefficient is the same as for the G V = 0 case and is just affected by the inclusion of vector interactions through the shift µ →μ rest 0 . As a consequence, the LP moves only along the chemical potential direction and the inhomogeneous phase is not suppressed [12].
Let us now investigate how the inclusion of vector interactions affects the competition between different kinds of modulations close to the LP. As already mentioned, throughout most of the inhomogeneous window the RKC at G V = 0 takes a sinusoidal shape. Indeed, performing a Fourier expansion of Eq. (12), one finds that the mass function is dominated by the lowest mode, unless the elliptic modulus is very close to ν = 1, which is only the case in the vicinity of the phase boundary to the homogeneous chirally broken phase [16]. In the following, we shall therefore consider a single cosine 1 as a prototype for real modulations with modulated density, namely we introduce M COS (z) = ∆ cos(qz) .
In order to investigate the competition between the cosine with corresponding modulated density and CDW solutions close to the LP, we start from Eq. (15) and plug into it the two different Ansätze, Eq. (19) and Eq. (11). After minimizing with respect to q in both cases, we can compute the GL free energy difference between the two different solutions. We find For G V = 0 this expression reduces to δΩ = (α 4 /4)∆ 4 [2]. Since α 4 < 0 in the inhomogeneous phase (otherwise nonzero mass gradients would not be favored), this shows that in the absence of vector interactions the cosine solution has a lower free energy than the CDW [2].
At the LP, on the other hand, α 4 vanishes. Very close to the LP δΩ is therefore dominated by the −R β β 42 term which is finite and positive for G V > 0 (cf. Eqs. (17) and (18)). We thus find that already for arbitrarily small values of G V there is always a region close to the LP where the CDW is favored over the cosine modulation.
Away from the LP, the other terms in Eq. (20) can play a role as well. The missing coefficients are Of course, we should be careful since the GL expansion breaks down as one moves away from the LP and the order parameters or their gradients become large. To get some idea about its validity, we also show the leading-order GL prediction for the phase boundary to the restored phase (dashed lines) in comparison with the exact numerical result. This seems to indicate that the boundary between the CDW and cosine regimes could be trusted in both examples. Note however that near the phase boundary the order parameter is always small and hence deviations can only be due to large gradients. In contrast, further away from the phase boundary the order parameter itself can be large as well. Therefore we refrain from making any quantitative statements in that regime.
In summary, the GL expansion suggests that, at least in proximity of the LP, the inclusion of vector interactions disfavors solutions with spatially modulated quark number density, making the CDW the thermodynamically favored solution.
It might be possible to extend the regime of validity of the present GL analysis further away from the Lifshitz point by following the procedure developed in [18]. A consistent In the following we focus again on the cosine ansatz as an approximation to the RKC and prototypical competitor to the CDW. In contrast to the latter, for which the density is spatially constant and therefore the ADA is exact, we expect the density of the cosine modulation to be inhomogeneous. Indeed, inserting the mass function Eq. (19) into Eq. (14) we find that in a leading-order (O(M 2 )) GL expansion δμ behaves like cos 2 (qz) ∝ (1 + cos(2qz)). For the shifted chemical potential we therefore consider the ansatz with a constant termμ 0 and a spatially modulated part with wave number 2q and amplitudẽ µ 1 . Accordingly, using Eq. (8), the density is given by with the spatial average n = µ−μ 0 2G V and an oscillating part with amplitude n 1 = −μ 1 2G V . Note that this behavior is qualitatively consistent with the analytically known density profile of the RKC at G V = 0 [12], which also varies periodically with half of the wavelength of the mass function.
Within this setup, the Dirac Hamiltonian H, Eq. (5), is found to be a non-diagonal matrix in momentum space, since momenta which differ by ±qe z or ±2qe z are coupled through the oscillating components of M andμ, respectively. The numerical determination of the eigenvalues is therefore rather involved but, exploiting the block structure due to the discreteness of the sets of coupled momenta, it is still feasible (for details see e.g. Refs. [1, 3,20]).
The resulting kinetic contribution to the thermodynamic potential, Eq. (4), is strongly divergent and needs to be regularized. Since sharp momentum cutoffs are known to fail when describing inhomogeneous phases [2, 20], we choose again a PV regularization scheme.
In order to stay as close as possible to the ADA calculations of Ref. [12], where the PV regularization was only applied to the G V -independent part H 0 after isolatingμ (see Eq. (9)), we now isolate the constant partμ 0 of the shifted chemical potential and apply the PV regularization only to the eigenvalues of the remaining part of the matrix. In other words, we write and then regularize Ω kin as whereẼ λ are the eigenvalues ofH, and sgn() is the sign function.
For this ansatz the thermodynamically favored state for the system is then obtained by solving for the stationary points of the thermodynamic potential 2 with respect toμ 0 andμ 1 , and minimizing it with respect to ∆ and q.
Our results for the order parameters and the density of the CDW and the sinusoidal modulation both within the ADA as well as with a modulated vector mean field are shown in Fig. 2 for G V = G S /2. For the cosine we find that the spatial average n of the modulated vector mean field agrees quite well with the density obtained in ADA. However, at the onset of the inhomogeneous phase, the amplitude n 1 of the density oscillations is almost as large as n , demonstrating that the ADA is a rather poor approximation in this region. Indeed, comparing the mass amplitudes and wave numbers of the ADA with those in the improved approach we find large differences, which go away only near the chiral restoration transition where the density oscillations become small.
We note that within the ADA the amplitude of the chiral condensate for the cosine ansatz jumps to a higher value at the phase transition from the homogeneous broken into the inhomogeneous phase. This feature disappears when we relax the ADA and consider a modulated vector condensate. In this case also the onset of the inhomogeneous phase is moved towards higher chemical potentials. By comparing the critical chemical potential with that of the CDW, it follows immediately that the CDW will become the favored modulation, at least in a certain region of the phase diagram. This is confirmed by comparing the free energies associated with the different modulations, as shown in Fig. 3  remain favored over the cosine all the way up to the chiral restoration transition. Whether this is the case is hard to tell due to our insufficient numerical resolution.

V. CONCLUSION
In this work, we investigated the effects of vector interactions on inhomogeneous chiral symmetry breaking within the NJL model, expanding on previous studies where the vector condensate was approximated by its spatial average [12]. In the present work we allowed for a more self-consistent treatment of the problem by considering spatially modulated vector condensates. Specifically, we studied a sinusoidal modulation of the order parameter in competition with a CDW, corresponding to a single plane wave. Aside from the known effect that vector interactions enlarge the inhomogeneous window in the model phase diagram [12], our main finding is that the sinusoidal modulation, which is energetically favored over the CDW if no vector interactions are present, eventually becomes disfavored when the vector coupling G V is increased. In particular, a Ginzburg-Landau analysis reveals that already for arbitrarily small values of G V , the CDW solution becomes favored in part of the inhomogeneous island in the model phase diagram. This result, which is expected to be valid in the proximity of the LP, is corroborated by a numerical study at zero temperature. This behavior can be understood by recalling that the repulsive vector interaction creates an additional energy cost if matter is not evenly spread out, so that any inhomogeneous distribution in the quark number density will become increasingly disfavored as the vector coupling becomes larger. As a consequence of this, the CDW solution, where the quark number density is homogeneous, becomes competitive with the originally more favored sinusoidal modulation, for which the density is spatially modulated.
The window where the CDW is favored over the sinusoidal modulation increases with G V .
This suggests that for sufficiently large vector couplings the CDW could be the most favored solution in the entire inhomogeneous phase. Whether this is indeed the case is however difficult to tell from our analysis since the numerical-diagonalization method suffers from a limited resolution, in particular near the phase boundary to the restored phase, while the GL expansion is valid only near the LP. The improved GL expansion developed in [18] might provide the tool required for this type of investigation.
We also note that, despite the improved treatment of the vector condensate, the sinusoidal ansatz we studied is not a selfconsistent solution. While at G V = 0 this ansatz is a very good approximation to the RKC, which in turn is a selfconsistent solution, we neither know the corresponding selfconsistent solution at G V = 0, nor whether it can still be well approximated by a cosine. It is therefore not excluded that with increasing G V the RKC could evolve continuously into a solution which persists to be energetically favored against both a sinusoidal shape and the CDW. It is also possible that there is a smooth transition between the (modified) RKC and a CDW at nonzero G V . In order to get more insight into these possibilities, it would be interesting to study more sophisticated ansätze with more independent Fourier modes for both mass functions and vector condensates.