{'fraction_non_alphanumeric': 0.04704608851326034, 'fraction_numerical': 0.014321995087222615, 'mean_word_length': 4.238822246455834, 'pattern_counts': {'":': 0, '<': 4, '<?xml version=': 0, '>': 7, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Multi Phase in Cold Dense Quark Matter

In this article, we study dynamic chiral symmetry breaking at zero temperature, finite chemical potential and external magnetic field with massless NJL model. We have proposed a mathematical method to classify phases in phase diagram of cold dense quark matter, and use mathematical analysis to identify the multi phase phenomenon among solutions for gap equation which means with fixed chemical potential and magnetic field, there could be two phases coexisting.


I. INTRODUCTION
The phase structure of QCD matter has always been an important and attractive topic in theoretical physics [1][2][3][4][5][6][7]. In relativistic heavy-ion collisions, the produced QCD matter will go though a phase transition or a crossover as time goes by. Either way, the state of QCD matter is believed to change from quark-gluon plasma to hadronic matter in this process. Its physical properties and dynamical behaviors such as chiral symmetry and confinement are altered along with the change of the state.
At the early stage of noncentral collision, the QCD matter produces extremely strong magnetic field [8][9][10][11], which brings about obvious magnetic effects. Moreover, the compact stellar objects such as magnetars are believed had strong magnetic field around 10 15 G at their surface [12,13], Therefore studying QCD matter's properties under the influence of magnetic field becomes a meaningful and important subject. So far, many relevant theories and models have been proposed and it is shown that the quark condensate are strengthened by magnetic field, which is known as 'Magnetic Catalysis' [14][15][16][17][18]. Consequently, the QCD phase diagram is related to magnetic field [19,20].
In this article, we will study the phase property of so called "cold dense quark matter", unlike the quark matter in high energy experiments, this matter has low temperature, hence we can establish models at zero temperature limitation. This kind of research could facilitate the study of compact stellar objects. Although similar projects have been thoroughly studied in early articles [21][22][23][24][25][26][27][28][29][30][31][32][33], but in this article we try to study phase diagram of cold dense matter in a different point of view, we have developed a mathematical method to classify phases in phase diagram, and find out new property from gap equation.
The model we employ is the two-flavor NJL model at chiral limitation with mean-field approximation [34][35][36][37], it is a good tool to investigate dynamical chiral symmetry breaking of nonperturbative QCD matter. One thing need mentioning here, in Asakawa's work [38], it was pointed that with the presence of chemical potential, the self-energy does not simply equal dynamical mass, which reveals with the help of the Fierz transformation. The actual self-energy should be written as Σ = σ + aγ 0 to guarantee the self-consistency of gap equation. In our case, this problem becomes much more complicate, the external magnetic field and chemical potential render self energy has four kinds of mean fields, Σ = σ + aγ 0 + bγ 5 γ 3 + cσ 12 . At zero temperature limitation, this will cause very chaotic situation in gap equation. But fortunately, except σ, the other three quantities are very small comparing to nonzero solutions of σ, we can ignore them for a schematic view of phase diagram.
The article is arranged as below, Section II is the deduction of gap equation with our new developed method different from Schwinger's "proper time", this method could handle more complicate models or ansatz such as Σ = σ+aγ 0 +bγ 5 γ 3 +cσ 12 , one can refer to Appendix A and B for more details of this method. Section III is the classification of phase diagram by using mathematical analysis and numerical analysis, in this section we mathematically define several areas in phase diagram which have different phase properties. Section IV is the conclusion.

II. GAP EQUATION
We start from the lagrangian below (which is the NJL model with mean field approximation) here A µ gives us the external magnetic field with the strength B that parallels to x 3 axis. q u , q d separately represent electric charge numbers of up quark and down quark, hence q u = 2 3 , q d = − 1 3 . We employ q f to generally represent q u , q d in following discussions, and the index 'f' can be assigned to 'u' or 'd'.
In the Lagrangian Eq. (1) we have employed the ansatz that Σ = σI 4 , while actually, with nonzero external magnetic field and chemical potential, Σ should have not only quark condensate but also vector, axial vector and tensor condensates. The last there kinds of condensates are thought to be too small to affect the general properties of QCD matter but some subtleties, therefore in this paper we assume Σ = σI 4 . Now we manage to get the fermion propagator with magnetic field and chemical potential, in the mean time, to identify the appropriate form of ε term (because it is zero temperature). According to Appendix A, in order to acquire the correct fermion propagator simply, one has to multiply Hamiltonian density with a factor (1 − iη), Consequently, the Lagrangian has changed to, the f here means we have separated flavor space, and the fermion field operator ψ in Eq. (4) is 4 components single flavor spinor rather than 8 components two flavor spinor in Eq. (1). Through the definition of partition function we have the effective action W ′ The principle of gap equation is to identify the least or local minimum value for effective action, which is hence we have in following discussion we can safely set factor (1 − iη) to be 1 at the left hand side of Eq. (9) with the limitation η → 0 + . Now we need to take care of TrŜ ′ f in Eq. (9), because there is chemical potential in the propagator, it is not convenient to use Schwinger's 'proper time' method. We have developed a new method to deal such situation. One can refer to our previous work [40] for a detailed introduction or to Appendix B for an overview.
According to Eq. (B5), we have Noticing, the complex factor in B ′ is (1 − iη) 2 rather than (1 − iη), because 2n|q f |eB ′ comes from the combination of quantization ofΠ 2 ⊥ and q f eBσ 12 which is the outcome of [Π f 1 ,Π f 2 ], both of them provide the factor (1 − iη) 2 . Secondly one should be aware due to nonzero chemical potential, Tr(Ŝ ′ f γ 0 ) is not 0, the consequence is that we should introduce a shift for chemical potential to prevent the inconsistency. But as we have mentioned at the beginning, such shift relates to a vector condensate and it is very small, we could exclude its effect in our qualitative results.
Continuing to adjust the expression of Eq. (10) Making a Wick rotation to Eq. (12), and applying proper time method Combining Eqs. (9) and (14), making a truncation to proper time 's' ( +∞ 0 ds → +∞ 1/Λ 2 ds), the gap equation could be simplified to (despite the trivial solution σ = 0) At the B → 0 + and µ → 0 + limits, Eq. (16) degenerates to the classic gap equation which we can use to determine the value of Λ and G [41], FIG. 1: In this figure, we treat σ as variable, eB and µ as preset parameters, and the vertical axis as function of σ.

A. The Boundary of Chemical Potential
There is a boundary to chemical potential. Before explain that, we have to vary gap equation Eq. (16) a little firstly The first row is the variation. In this way, we can study two functions f (σ, eB) and h(σ, µ, eB) separately beside gap equation. From the expression of h(σ, µ, eB) in Eq. (19) we know if chemical potential µ is smaller than dynamic mass σ, then the gap equation is simplified to there will be no µ-dependent dynamic mass in existence, we define these ordinary 'only-magnet-dependent' dynamic mass (ODM) as σ 0 . σ 0 roughly gives us lower boundary of µ that the solutions of gap equation Eq. (19) are beyond 'ordinary'. In the matter of fact, the actual lower boundary is generally smaller than σ 0 with all eBs, take Fig. (1) for example, we have three solutions for gap equation in the case of µ 1 , σ 0 is the ODM, σ ′ is also a valid solution for Eq. (19) while σ ′′ is not, because from mathematical analysis we know at σ ′′ the effective action Eq. (7) has a local maximum rather than minimum, it is not what we need. In following discussion we will ignore these maximum points.
The contacting points of f (σ, eB) and h(σ, µ, eB) in Fig. (1) is not a valid solution neither, but it is a boundary point for chemical potential. Comparing µ 1 line with µ 2 , we can see when µ > µ 2 , Eq. (19) has valid solutions beyond ODM, and analytically it is true, because h(σ, µ, eB) is a monotonically increasing function by µ, as long as µ is bigger than a specific value at different eBs, the solution will be not just ODM. Therefore the real lower boundary of µ, defined as µ low , is always lower than σ 0 . To identify µ low , we have to solve the simultaneous equations of σ, Beside lower boundary, we also have upper boundary for chemical potential, defined as µ up , when µ exceeds such boundary, Eq. (19) has no valid solution (f (µ, eB) has no intersection with h(σ, µ, eB), or the intersection point represents maximum rather minimum of effective action), which means chiral symmetry restores, Wigner phase arises.
For the function f (σ, eB), by L'Hôpital's rule, therefore we have here '∼' reads "f (σ, eB) is asymptotic to f 0 (σ, eB) as σ tends to 0 + " [42]. If someone get confused on the units problem in the logarithm function of σ, one can use to replace the one in Eq. (24). For the function h(σ, µ, eB), Noticing, as σ → 0 + , both f (σ, eB) and h(σ, µ, eB) have the same dominant asymptotic behaviors, closing to (−2eB ln σ), this renders quite complicate relations between them, and the solutions too, when they are close enough. But still, we can make a use of their asymptotic behaviors to 'roughly' identify the upper boundary of chemical potential by equation of µ, the ln σ terms in Eq. (26) are perfectly canceled on both side, which leaves us an implicit function with respect to µ up and eB.
In the above discussion we have mentioned "roughly identify the upper boundary", that's because the actual upper boundary is beyond Eq. (26). There are two kinds of situations. First, f (σ, eB) is a smooth function of σ but h(σ, µ, eB) is not, beside that, f (σ, eB) and h(σ, µ up , eB) are not only very close as σ → 0 + , but also close enough when σ stretches to a finite value, say 0.1GeV, these two reasons cause multiple intersections before they distinctly separate, and of cause some of the intersections account for nonzero dynamic mass, e.g. Fig. (2). But that doesn't say Eq. (26) is of no use, we can add a tiny modification to µ up from Eq. (26), which brings us the actual upper boundary µ ub = µ up + △µ, in the case of Fig. (2), △µ ≈ 0.0018GeV, it is a small quantity, therefore we can say Eq. (26) roughly identifies upper boundary.
The second situation, through numerical result, we find that µ up decreases along with increasing magnetic field, while the ODM σ 0 is increasing, therefore when µ up , as the function of eB, exceeds a threshold, f (σ, eB) will always σ0 is a valid intersection between f (σ, eB) and h(σ, µup, eB), but is not valid between f (σ, eB) and h(σ, σ0, eB).
intersect with h(σ, µ up , eB) at a ODM point, e.g. Fig. (3), this implies that the actual upper boundary, defined as µ ub , depends on σ 0 rather than µ up .
So far, we have discussed the properties of σ 0 , µ low and µ up , we put them in Fig. (4). The modifications of µ up is also sketchily plotted in this figure as error bars, a precise demonstration of these modifications is shown in Fig.  (5), which, as we can see, is relatively small and ruleless, therefore for a qualitative discussion, we can just talk about µ up instead of µ ub when µ ub = µ up + △µ. Noticeably, when magnetic field is not strong enough, relatively, the upper boundary of chemical potential, which separates chiral restored phase and chiral broken phase, is highly ruleless, it has irregular oscillation, but when eB keeps increasing, µ up begins regularly decreasing, we will discuss this in the conclusion.
From Fig. (4), the three functions of eB have roughly divided the diagram into several areas, this gives us phase diagram with nonzero chemical potential and magnetic field at zero temperature, as shown in Fig. (6). The diagram is mainly divided to three kinds of areas, the ordinary Nambu phase, the multi phase area and the Wigner phase, and in multi phase area (MPA), it is divided into three subareas by function min(µ up , σ 0 ). The ordinary Nambu phase has only one solution, the ODM, determined by Eq. (20), it depends only on magnetic field, has no relationship with chemical potential. The Wigner phase is the chiral restored phase, it always has zero dynamic mass σ = 0 in this model (chiral limit NJL model) we study. The most interesting part is the MPA, its three subareas have different properties of the solutions of gap equation Eq. (19). We will discuss this in detail next subsection.

B. The Solutions in Multi Phase Area
In order to clarify the phase properties in MPA of Fig. (6), we have to treat h(σ, µ, eB) in a proper way. Due to the Heaviside step function, h(σ, µ, eB), as a function of σ, is continue but not smooth, it has stages, and at the lowest stage, h becomes a constant 4π 2 /G, therefore we consider h(σ, µ, eB) as a two sections function, section I, the stairs (h > 4π 2 /G), in this section, h could have many stairs, inside the stair, h is smooth, at the point of two FIG. 4: σ0, µup, µ low are implicit functions of eB. The error bars of µup represent modification △µ to µup, they are only sketch, the actual modifications are relatively small quantities. When eB > 0.16GeV 2 , µup < σ0, the actual upper boundary of chemical potential depends on σ0 rather than µup. Along with eB's increasing, µ low and µup are decreasing, and getting closer and closer, this property could be easily understood through mathematics analysis. Noticing, when eB ranges from 0 to 0.15GeV 2 , µup has irregular oscillation.
FIG. 5: The modifications of µup, lead to actual upper boundaries of chemical potential. The range of eB is from 0 to 0.16GeV 2 , because when eB > 0.16GeV 2 , the actual upper boundaries is determined by σ0, there is no need to consider modification for µup. The biggest modification happens at eB = 0.07GeV 2 , △µ ≈ 0.017GeV, while the smallest modification, except 0, is △µ ≈ 0.0003GeV at eB = 0.03GeV 2 . At some points, △µ = 0, no modification is needed, µ ub = µup. stairs contacting, h is continue but not smooth, the h function in Fig. (2) is a good example, section II, the ground (h = 4π 2 /G), of cause at the point of stairs and ground contacting, h is continue but not smooth, too. When f and h intersect at the ground, we have ODM, when at the stairs, we have UDM (unordinary dynamic mass).
MPA I is a area that only has UDMs. It quite clear that in this area no ODM solution is allowed, because µ > σ 0 . In MPA I, we also have µ > µ low , µ low means the point of first contact of f and h, when µ exceeds the first contact point, f and h will have at least one intersection, of cause an UDM, or multi UDMs simultaneously, e.g. Fig. (7). But multi UDMs does not happens all the time, generally, when the solution of gap equation is closing to a contacting point of two stairs of h, we can find another valid solution in the other stair. The σ-µ relation could be like Fig. (8), which is just a small part of a big picture. When µ ranges from µ low to µ up , the cascade could happen many times, and the multi phase areas only exist for a short range.
MPA II is the area that definitely has multi phase simultaneously, in this area, we always have at least two valid solutions from the gap equation, and one of which is an ODM. Because µ < σ 0 , f has an intersection with h at the ground, while µ low < µ < µ up , f always cuts through the stairs of h, that brings us another valid solution. The f and h's relation is quite like Fig. (7), except the lowest stair of h is ground 4π 2 /G. The σ-µ relation about this area is shown in Fig. (9).
MPA III is just another ordinary Nambu phase area, in this area, all dynamic mass are determined by Eq. (20). From Fig. (6), it seems MPA III and ordinary Nambu phase from the bottom are separated by MPA II, but considering µ up and µ low are moving really closer with increasing eB, maybe they will connect when eB is strong enough.

IV. THE CONCLUSION
In this article, we have developed a mathematical analysis method to help drawing phase diagram of cold dense quark matter, and roughly divide phase diagram into several areas, each area has a unique phase property. The phase diagram and division depend on three kinds of quantities, σ 0 , µ low , µ up , they are all implicit functions of eB, the equations that depict these quantities are separately Eqs. (20), (21), (26). σ 0 is the ordinary dynamic mass in NJL model with external magnetic field, when chemical potential is involved, there are some dynamic mass deviating ODM, therefore σ 0 is the base line of all other phases. µ low is the first contacting point of f (σ, eB) and h(σ, µ, eB), when chemical potential is smaller than µ low , f and h can only have the ODM solutions for gap equation, therefore µ low is the dividing line of ordinary dynamic mass and unordinary dynamic mass. The meaning of µ up is a bit complicate, generally speaking, it is the dividing line of multi phase area and simple phase area, when magnetic field is relatively weak, say 0.15GeV 2 comparing to 0.3GeV 2 , µ up separates the multi phase area from Wigner phase area, when magnetic field is strong, µ up separates multi phase area from another ordinary Nambu phase area. Of cause this classification is not as subtle as the phase diagram in [33], but all these dividing lines are mathematically definable.
When magnetic field is below about 0.14GeV 2 , the upper boundary µ up of chemical potential has an interesting irregular oscillation, this could be caused by quantum fluctuation around the critical point. However, a strong enough magnetic field would smear the fluctuation, that's why when eB > 0.14GeV 2 , µ up has a regular smooth descending. Mathematically speaking, when eB is small, there are several Landau levels play roles in the gap equation due to Heaviside step functions in h(σ, µ, eB), these levels cause the irregular solutions in Eq. (26), but when eB is strong enough, that leaves us only the lowest Landau level (θ(µ − σ) term in h), and it is much more regular and predictable. We can roughly estimate from which point on µ up becomes regular, firstly, the phase transition line between chiral restored phase and chiral broken phase is around 0.3GeV, so we assume µ = 0.3GeV, and now if we want only lowest Landau level involving in h, it requires (µ − 2|q d |eB) ≤ 0, which leads to eB ≥ 0.135GeV 2 , it is pretty close to 0.14GeV 2 .
In this article we have proved the existence of multi phase, shown in Figs. (8) and (9). The case in Fig. (8) belongs to MPA I of Fig. (6), in this area, not any chemical potential guarantees two valid solutions to the gap equation, only when chemical potential belongs to some specific intervals, the two valid solutions can be found. When we use numerical method to solve Eq. (19), if one solution is close to the joint point of two different stairs in h(σ, µ, eB) such as µ 2 − 2n|q f |eB, then maybe we could find another solution on the other side of the joint point. It probably has no mathematical equations to identify whether or not there is multi solution around σ = µ 2 − 2n|q f |eB, we suspect this phenomenon happens whenever dynamic mass jumps from one stair to another with chemical potential changing, and the smaller dynamic mass is, the shorter intervals chemical potential need going through. The case in Fig. (9) belongs to MPA II, as we can see in Fig. (6), this area spread from weak magnetic field to strong magnetic field, which means multi phase formed by ODM and UDM always exists as long as magnetic field is nonzero, this is understandable, because no matter how strong or weak the magnetic field is, from Nambu phase to Wigner phase, all particles have to pass through lowest Landau level, and this level causes multi phase.
The physical effect of multi phase is energy level transition, we take the cases from Figs. (8) and (9) as examples, fix magnetic field and chemical potential at specific values which guarantee multi phase, treat σ as free variable of free energy density F (W = F d 4 x) from Eq. (7), unsurprisingly there are two minimums of F in the intervals we choose, seen in Fig. (10). Because σ 1,2 (or σ 0,3 exist in the same external conditions (chemical potential and magnetic field), they can transfer to each other accompanied by energy absorption or radiation. If the absorption happens, absorbing photon for instance, some particles jump to higher Landau level (only one level higher for sure, because multi phase happens between adjacent levels), or vice versa. Interestingly when energy level transition happens, not only Landau level, but also dynamic mass changes, phenomenally speaking, the "structure" of quark matter has changed. We know if a thermal system is in equilibrium state with multiple phases, all phases must fulfill three kinds of equilibrium, thermal equilibrium (equal temperature), mechanical equilibrium (equal pressure) and diffusive equilibrium (equal chemical potential). In the multi phase case here, the quark matter remains at zero temperature (thermal equilibrium), with fixing chemical potential (diffusive equilibrium) and magnetic field (it is an external condition and space-timeindependent), the free energy density (or pressure, P = −F ) is not equal, mechanical equilibrium is not fulfilled, therefore the particles in high energy state such as σ 1 in Fig. (10) will be pushed away by or transfer to particles in lower energy state σ 2 , the multi phase state could not stably exist.
In this article, we employ zero temperature limitation and the ansatz Σ = σ, but as we mentioned at the beginning, a complete self energy should involve four kinds of condensate, if considering all these condensates, undoubtedly properties of phases for gap equations should be more complicate and abundant, in following work, we would like to do more detailed study, but as to a schematic view of quark matter at zero temperature, the ansatz in this article is adequate. The advantage of zero temperature limitation is that gap equation can be more "clear", and we are able to rely on mathematical tools to analyse gap equation, but when temperature goes to nonzero, we have to depend on numerical methods, it could be hard to find some properties of the solutions such as multi phase. We believe when temperature is low, multi phase phenomenon could still exist, because from zero temperature to low temperature, it is a continuous process, the properties of phases should be continuously varying. Also in following studies, we would like explore the cases at high temperature, although the particles' thermal motions at high temperature will smear many effects that can be found at low temperature, but if the multi phase is still there, that could be significant.
To demonstrate this method we have developed is valid, here we consider another example, the zero temperature and finite chemical potential case of free fermion propagator.
In momentum space, we have this looks a bit different from what we need. First of all, in the momentum integral, we can always make a shift to p 0 , p 0 + µ → p 0 . Secondly, in the denominator of Eq. (A9), p 2 0 − ω 2 implies the poles which require p 0 = ±ω, we could use this relation to replace ω in ε term, which gives us .

(A10)
That is what we need. Of cause in ε term, the replacement of ω with p 0 seems a little undemanding. One can also prove rigorously from complex analysis that such replacement is legitimate.