Anatomy of the dense QCD matter from canonical sectors

We investigate the nuclear and the quark matter at finite real chemical potential ($\mu_\mathrm{R}$) and low temperature from the viewpoint of the canonical sectors constructed via the imaginary chemical potential region. Based on the large $N_\mathrm{c}$ estimation, where $N_\mathrm{c}$ is the number of color, we can discuss the confinement-deconfinement nature at finite $\mu_\mathrm{R}$ from the canonical sectors. We found the expectation that the sharp change of canonical sectors at $\mu_\mathrm{R} \sim M_\mathrm{B}/N_\mathrm{c}$, where $M_\mathrm{B}$ is the lowest baryon mass, is happen in the large $N_\mathrm{c}$ regime, and it is matched with the quarkyonic picture. In addition, we discussed the color superconductivity and the chiral properties from the structure of canonical sectors. Even in the present anatomy from the canonical sectors, we can have the suitable picture for the dense QCD matter.


I. INTRODUCTION
Exploring the phase structure of Quantum Chromodynamics (QCD) at finite temperature (T ) and real chemical potential (µ R ) is an important and interesting subject in the several research fields such as elementary particle, hadron and nuclear physics. Particularly, the moderate µ R region attracts much more attention recently because some exotic phases such as the color-superconducting, quarkyonic and inhomogeneous chiral symmetry broken phases are expected to be appeared in the region; for example, see Ref. [1]. If we can obtain the QCD phase diagram starting from the first principle calculation that is the lattice QCD simulation, there are no unclearness, but it is not feasible at moderate µ R because of the well known sign problem; see Ref. [2] for a review of the sign problem and Refs. [3][4][5][6][7][8] for recent progresses in methods to tackle the sign problem as an example. Therefore, several expectations at moderate µ R have been obtained by using QCD effective models; see Refs. [1,9] as an example.
In the early stage of the study for the QCD phase diagram, the first-order transition was expected to be appeared at moderate µ R even with sufficiently low T . However, the duality between the hadron phase and the color superconducting phase, which include not only the color-flavor locking (CFL) phase but also the twoflavor color-superconducting (2SC) phase, has been proposed [10][11][12]; there are one to one correspondence of elementary excitation modes between phases. It is the so called the quark-hadron continuity. In addition, there are several studies that predict the crossover between the hadron phase and the deconfined quark matter by using the QCD effective model [13] and also the Ginzburg-Landau analysis [14]. Of course, the duality is not confirmed yet; see Refs. [15] for recent progresses.
To understand the QCD phase diagram at moderate µ R , the quakyonic phase [16] may play a crucial role. * kashiwa@fit.ac.jp † kounoh@cc.saga-u.ac.jp The quarkyonic phase is first proposed in the large N c QCD by using the N c counting where N c is the number of color; see Refs. [17,18]. In the large N c , quarks can be treated as the probe if µ R is not reached to O(N c ). Then the pressure is O(N 0 c ) in the confined phase because the physical degree of freedoms are the gluballs and baryons, but it does O(N 2 c ) in the deconfined phase because physical degree of freedoms are gluons and quarks. Interestingly, the quark number density start to have nonzero value when µ R reaches to M B /N c where M B is the lowest baryon mass. In this case, the pressure is O(N 1 c ), and it can be interpreted as follows: The thermodynamic quantities are dominated by quarks inside the Fermi sea, but the physical excitation modes on the Fermi surface are corresponding to baryonic because the confined nature is not changed; the quarks are probe and thus they can not modify the gauge field configuration. Recently, the quarkyonic phase has also been investigated by using the top-down approach based on the AdS/CFT correspondence [19]. Unfortunately, the quakyonic phase is not clear when we consider small N c such as N c = 3, but some discussions have been done with QCD effective models; for example, see Ref. [20]. There is the discussion that the quarkyonic phase still important in the N c = 3 system such as the new confinement-deconfinement transition scenario [21]. In addition, the quarkyonic phase may affect the neutron star properties such as the massradius relation [22].
To understand the physical degree of freedom, the canonical ensembles may provide important information since they are directly related with each quark number. Therefore, we employ the canonical ensemble method [23][24][25][26][27] to discuss the QCD phase structure at moderate µ R with low T in this paper. The canonical sectors are constructed by using the imaginary chemical potential (µ I ) and then several knowledge of QCD at finite µ I play a crucial role; see Refs. [28][29][30][31][32][33][34] as an example. In this paper, we first consider the quakyonic phase. The color superconducting phase is, of course, interesting phase on the QCD phase diagram, but the gauge symmetry issue in the canonical ensemble method is a tricky issue. Thus, it is difficult to investigate the phase in the present approach, but we show some qualitative discussions for the color superconducting phase with some ansatz in this paper. In addition, the chiral symmetry restoration with increasing µ R is discussed from the canonical sectors. This paper is organized as follows. In the next section II, we explain the formulation of the canonical ensemble method. Section III explains the discussions about the canonical sectors in the large N c and Sec. IV shows the simple estimation of the canonical sectors by using the Polyakov-loop extended Nambu-Jona-Lasinio model as an example. The chiral properties and the color superconducting are discussed in Sec. V and VI, respectively. Section VII is devoted to the summary.

II. CANONICAL ENSEMBLE METHOD
Throughout whole discussions in this paper, we consider large but finite size continues system because the thermodynamic limit requires careful treatment for the infinite sums. Starting from the grand-canonical partition function (Z GC ) at finite T and θ := µ I /T , we can construct the canonical partition function (Z C ) with fixed quark number (Q) as where H menas the Hamiltonian, k ∈ Z, z = exp(2πi/N c ) is the Z Nc factor andn is the quark number operator. In this paper, we do not explicitly show T for the argument of the partition function and also other quantities because we are interested in µ R effects. With the fugacity expansion, we have where n B = n/N c is the baryon number and µ B = N c µ R does the baryon chemical potential; we here use the fact that N c multiples of n only contribute Z C because of the Roberge-Weiss (RW) periodicity; see Eq. (1) and Ref. [28] as an example. By using the above relations and the lattice QCD data at finite µ I , we can investigate the QCD phase structure at finite µ R with certain T where the numerical error induced from the Fourier transformation can be controlled [23][24][25][26][27]. It is noted that above expression is valid even for low T and high µ R ; numerical confirmations of it in QCD effective models can be seen in Refs. [35,36]. Therefore, Z GC (µ R ) is constructed from Z GC (θ), and vice versa. Finally, we have the inverse relation as This relation means that Z GC (θ) is decomposed by the canonical sectors (oscillating modes) and thus the canonical sectors survey elementary excitation modes via the oscillating behaviors. This fact has been used to investigate the confinement-deconfinement nature at finite T with µ R = 0 [33].

III. LARGE Nc ESTIMATION
Via the N c counting of meson, baryon, gluball, gluon and quark contributions, the quakyonic matter has been proposed in the large N c QCD; see Fig. 1. Since quark loops are 1/N c suppressed, the oscillation of quantities at finite θ below the critical temperature T c may be approximated by simple cos(N c θ) function. Therefore, we can assume where a and b 1 should depend on T and the spatial volume (V ). For example, the similar functional form of Eq. (4) is obtained in the strong coupling calculation; see Refs. [37,38]. With decreasing T even at small N c , the oscillation in terms of θ becomes weak and thus b 1 → 0 with T → 0 because θ can be transformed as the temporal boundary condition of quarks. Higher order terms are expected to be suppressed and thus we neglected here. The first term in Eq. (4), a, represents the gluball contributions because it does not have the µ-dependence. The second term should mainly contain the baryonic contributions.
In this work, we consider large but finite N c , 1 N c < ∞, from following reasons: The phase rotation of the Polyakov loop in terms of θ does not happen in the large N c limit because Z Nc broken quark contributions, which induces the RW periodicity, can not be modified Z Nc symmetric gluon contributions. This indicates that the period in terms of θ is 2π at least in the high T region as mentioned in Refs. [39][40][41], which is shape contrast with the finite N c case, when quarks become the probe. At low T , the system is perfectly dominated by gluballs and baryons in the large N c limit and then the periodicity . Solid lines represent the phase transition lines and closed circles mean the critical endpoint. The legend CSC means the color superconducting phase such as the CFL and 2SC and TRW denots the Roberge-Weiss endpoint temperature which is almost corresponding to the deconfinement critical temperature at µ = 0 with large Nc. In the bottom panel, the confined phase is not clear meaning above the liquid-gas transition which is not explicitly shown in the figure. Also, we assume that there is no chiral phase transition at low T in the figure; it is, of course, not confirmed yet.
issue may be relaxed; there is the RW periodicity. However, the period of the RW periodicity at low T becomes 0 in the large N c limit and then it is highly nontrivial that the Fourier transformation is well defined or not when we consider the order of operation about N c → ∞ and the integration of Fourier transformation; see Sec. II. Therefore, we need the infinitesimally small but nonzero back-reaction of quarks to the gluon contributions in this work.
Via the change of variables in the partition function, θ can be absorbed into the temporal boundary condition of quarks and thus it is irrelevant at low T . This indicates the important consequence: If the system is confined at µ = 0 with fixed T , the imaginary chemical potential re-gion should be the confined region. Then, the canonical sectors constructed by using the imaginary chemical potential should be independent of µ R and thus all canonical sectors only have the confined information it should be. This fact is valid even if N c is small if the confined nature is strong enough such as the T ∼ 0 situation. Below, we consider sufficiently low T because both the fugacity and Z C (k) with k = 1, 2, · · · depend on T and thus discussions about the T -dependence is very difficult unlike the µ R -dependence.
The canonical sectors with Eq. (4) become because . (6) These results are, of course, consequences from the properties of Fourier transformation; we put them here to show the consistency. Then,we have Below, we consider large N c and nonzero positive µ R and thus we neglect the third term. Since the second term depends on µ R and thus there is the region where the first and the second term are balanced atμ R ; Since the second term is the N c quarks (baryon) contribution and thus this energy scale may be related with the quarkyonic phase transition. It is well known that the quakyonic phase transition is happen around When we match this value with Eq. (8), we have Therefore, we obtain This result seems to be natural because of following two reasons. First, b 1 comes from the baryon contribution in the case and thus appearance of M B is natural. Second, the oscillation of Z GC in terms of θ should be vanished with T → 0, and this fact is realized in Eq. (10) as b 1 → 0 with T → 0. In addition, we can expect this kind of the functional form is qualitatively obtained in the QCD effective model with suitable simplifications; see the next section IV.
The pressure (p) is obtained from the thermodynamic relation as where β is the inverse temperature, β = 1/T . The b 1 term induces the N 1 c -order contributions because M B ∝ N c until µ R reaches N 1 c -order. The Z C (0) term does not have the µ R -dependence but e NcµR/T Z C (N c ) have the dependence and thus the b 1 term leads the opening contribution to nonzero quark number density (n q ); Therefore, n q ∼ 0 remains until µ R ∼ M B /N c and it is turned into N 1 c -order above the value with large N c and low T . This behavior is also consistent with the quarkyonic picture in large N c QCD about the quark number density.
Finally, in this section, we discuss following two scenarios. Via the present treatment used in this study, at least, we discuss precursory phenomenon for the phase transitions. The first one is (A): the change of infinite tower of the canonical sectors with increasing µ R and the other is (B): the nonexistence of the change.
. Since a should be the gluball contribution and thus it leads N 0 c contributions to the pressure. Interestingly, there is the balanced regionμ R − < µ c <μ R + with ∼ Λ QCD /N c where we use Λ QCD in the order counting. If each coefficient has, at least, the 1/n suppression factor, the dominant canonical sectors are changed as where the infinite tower of the canonical sectors are spanned with increasing µ R because of the exponential suppression factors for each canonical sector. Since quark loops are suppressed in the N c → ∞ limit and thus we can naturally expect such 1/n suppression factor; we may expect stronger suppression factor from the lattice QCD data at finite T [42]. If we have stronger suppression factor, discussions in this paper are unchanged.
In the large N c , the balanced region is very tiny, but it is enlarged when N c becomes smaller. This indicates that there is no clear phase transition in the realistic N c at moderate µ R with low T . This expectation is consistent with recent discussion for the moderate µ R region; there is no phase transition between the nuclear and quark matters. In contrast, in the case of realistic N c case, it is difficult to obtain clear consideration, but we can expect the balanced region is enlarged from the above estimation of the region and thus we do not have clear quakyonic phase transition energy scale; see the bottom panel of Fig. 1 as a schematic figure of the QCD phase diagram. This scenario is most likely scenario in QCD.

Scenario (B):
It should be noted that all higher-order canonical sectors start to contribute the system at µ R if each higher-order canonical sector does not have the 1/n-type suppression factor. Then, we can not have the infinite tower picture unlike the scenario (A). In this case, we may have the sharp energy scale of the quarkyonic phase transition even in the realistic N c = 3 system form the comparison between Z GC (0) and Z GC (N c ). However, we can find the suppression factor from the simple estimation by using the QCD effective model as in the next section IV and thus this scenario is not likely scenario in QCD. It should be noted that Z GC is not simply converged in this scenario with increasing n and thus it seems to be rejected from the mathematical sense.

IV. SIMPLE MODEL ESTIMATION
In this section, we show canonical sectors estimated by employing the Polyakov-loop extended Nambu-Jona-Lasinio (PNJL) model [43,44] as an example. The PNJL model is the extended model of the NJL model to include the Polyakov-loop dynamics by considering the meanfield of the temporal component of the gluon field (A 4 ) with the homogeneous ansatz. Of course, it does not contain the exact confinement mechanism, but it can mimic not only the approximated confinement nature but also several important properties such as the RW periodicity and its transition which are closely related with physical degree of freedom of the system [45,46]. Therefore, we here employ the PNJL model to estimate canonical sectors and show the estimation is matched with the result obtained in the previous section III.
In the following estimation, we assume that the constituent quark mass (M ) is order of the QCD energy scale (Λ QCD ) in the finite θ region with low T when we estimate the Fourier transformation. In the PNJL model, we usually employ the Polyakov-gauge fixing, ∂A 4 = 0, and then A 4 is diagonalized by using the gauge degree of freedom remained in the spatial component. The Lagrangian density of the two-flavor PNJL model is given by where m 0 denotes the current quark mass, D µ is the covariant derivative D µ = ∂ µ + iδ 4 µ A µ , G is the coupling constant and U does the gluonic contribution. In this paper, U is not important and thus we do not show the explicit functional form; for example, see Refs. [43,47,48] With the mean-field approximation, the effective poten-tial is given by where N f = 2,μ = µ + iA 4 and here E ∓ = E ∓ µ, σ = qq and E = p 2 + M 2 with M = m − 2Gσ. The definition of the Polyakov loop (Φ) and its conjugate (Φ) in the model are where tr c is the trace acts on the color space. The above effective potential (16) is corresponding to the leadingorder result of the 1/N c expansion and then the higher order contributions which are corresponding to the meson loop contributions are neglected; see Ref. [49] for the result in the NJL model as an example.
Since we here consider the low T region and thus the Polyakov-loop dynamics is usually decoupled from the system; we here consider the T region where the Polyakov loop and some other higher-order loops are sufficiently weaker than exp(−N c β(E ∓ µ)). Therefore, the effective potential (15) may be simplified as wheref It should be noted that we here consider finite imaginary chemical potential and thus µ = iµ I . When we consider N c > 3, there should be contributions of higher order loops [50] and thus f ∓ must be modified comparing with those in N c = 3, but we can expect the appearance of the e −Ncβ(E+µ) term in the logarithmic part because of properties of SU (N c ) Lie group; it is because ln det term in Eq. (15) must have exp[−N c β(E ∓ µ)]. This term mimics the N c -quark state and the suppression of other m-quark states with m = 1, · · · , N c − 1 at low T does the QCD confinement nature; see Ref. [51] for discussions in the case of N c = 3. In addition, T is set to be sufficiently small and thus we can approximate Eq. (15) as where σ depends on θ, but the dominant part can be expected as θ-independent and thus we approximate the second and third terms in Eq. (18) into the θ-independent a. In the case of the large N c , N c M dominates the integration (static limit) as because the exponential suppression is much stronger than p 2 in the large N c , and where and˜ are the infinitesimal value proportional to 1/N c and the infinitesimal value which manifests˜ 1/N c . Therefore, we consider the momentum integration up to the region where we do not reach¯ scale because of the strong exponential suppression effect. Thus, we have where a is the θ-independent and b does the θ-dependent coefficients. Here, we use M 2 /M 2 + 1 ∼ M . If the saddle-point approximation is good, we can estimate the grand-canonical partition function as Therefore, we finally obtain Z ∼ e −βV a 1 − βV b e −NcM/T cos(N c θ) + · · · =ã +b e −NcM/T cos(N c θ) + · · · , with sufficiently low T and fixed V . The second term is corresponding to b 1 in the previous section. Interestingly, although we use too much simplification, we can find the exp(−N c M/T ) ∼ exp(−M B /T ) factor in the present estimation. This form is exactly the same one that we obtained in the large N c estimation via the matching with quarkyonic phase transition energy scale in Sec. III; the difference of coefficients between the present and previous estimation can be adjusted by µ c → µ c + α and α may be order O( ln Nc Nc ). In principle, higher-order contributions can exist in Eq. (25) and then we have exp(−mN c M/T ) terms with m = 1, 2, · · · ; those are not important in the region where the dominant oscillating mode is cos(N c θ) when µ R is not large enough. It should be noted that we expand the logarithmic functions by exp(−N c M/T ) in Eq. (20) and thus higher-order expansion terms should have, at least, the 1/n suppression factor; existence of this factor is assumed in the previous section, but we can find it in the present estimation.

V. CHIRAL PROPERTIES
In this section, we discuss the chiral properties from the canonical sectors. In the case of the chiral condensate, it can have finite value at finite θ at low T . Then, the chiral condensate at finite µ R can be expressed as At low µ R , higher-order canonical sectors than the n = 0 sector are irrelevant and then the contribution of σ C (0) dominate σ (µ R ). It is known that the absolute value of the chiral condensate at finite θ becomes larger than that at µ = 0 from the lattice QCD simulation [29] and the QCD effective model calculations [52]. This indicates that the oscillating behavior at sufficiently low T is expected as where a σ and b σ are positive coefficients and σ means the absolute value of the chiral condensate. In the T → 0 limit, a σ dominates the system and then we can expect a σ is the Λ QCD -order; the first term σ C (0) can be assumed as the Λ QCD -order in the setup. When where sgn is the sign function. Because of the sign difference, the n = 0 canonical sector is weaken by the n = N c canonical sector. If µ R increases more and more, higherorder canonical sectors join the cancellation of the n = 0 canonical sector and then the sign for each canonical sectors become important. To discuss more details, we need actual numerical data and thus it is left unclear problem which will be clarify in our future work. Present discussion is the picture of the chiral symmetry restoration with increasing µ R at sufficiently low T from the viewpoint of the canonical sectors.
If the chiral phase transition which can be characterized by the corresponding local order-parameter happens at µ R = µ σ , Z GC becomes zero very close to µ σ on the complex µ plane; for example, this fact plays a crucial role in the Lee-Yang zero analysis [53]. With sufficiently low T , we can estimate the transition point from and then where M is the constituent quark mass which depends on the current quark mass and we here temporally introduce the infinitesimally small µ I = .
It should be noted that we here consider the sufficiently large but finite size system and thus phase transition is smeared. It is, however, Lee-Yang zeros can appear very close to the real axis and thus we can evaluate the phase transition point from the finite size system without any inconsistencies.
In the large N c with sufficiently low T , we can consider following three scenarios (I)-(III).

Scenario (I):
If the second term is 1/N c suppressed in Eq. (29) and the cos(N c θ) function dominates the oscillating behavior of Z GC (θ), the transition point is coincident with the value of the quarkyonic phase transition energy scale. This is corresponding to the situation that A and B 1 are the same order.

Scenario (II):
If ln(A/B 1 ) in the second term is proportional to N c , the transition point is shifted from M .

Scenario (III):
If ln(A/B 1 ) in the second term is proportional to N 2 c or higher, the transition point is diverged. At low T , this scenario seems as the unfeasible scenario because of the following discussion based on the physical degree of freedom.
Since the Stefan-Boltzmann limit of the pressure is proportional to N 2 c because of the gluon degree of freedom, but it should be lower number in the sufficiently low T region because the physical degree of freedoms are gluballs and baryons even at finite θ. Therefore, in the present energy scale, the scenario (I) or (II) is realized because the pressure and the partition function are related with each other via Eq. (11); for example, see Ref. [54] for the normalized pressure on the lattice by that with the Stefan-Boltzmann limit for several N c . These discussion, of course, must be modified with larger T .
Above discussions are implicitly assumed that condensates are homogeneous. Even if we consider the spatial inhomogeneity, above discussions are valid if we replace σ as σ(x). Then, we can discuss the inhomogensously chiral-symmetry broken phase which is the so called the real kink crystal [55] by introduce the external field which breaks the transnational invariance (J δ ) and take the J δ → 0 limit. In addition, via the same procedure as in the diquark condensate which will be explained in the next section, we can consider the condensation of the neutral pion π 0 (x) (µ R ) because the pion condensate does not appear in the finite θ region in the J δ → 0 limit. Then, we can consistently investigate the dual chiral density wave [56] which is another type of the inhomogeneous chiral-symmetry broken phase in addition to the real kink crystal.

VI. COLOR SUPERCONDUCTIVITY
Here, we discuss the color superconducting from the viewpoint of the canonical sectors. However, the color superconductivity is 1/N c suppressed [57] and thus it is difficult to obtain clear picture from the present canonical sector approach without actual numerical calculations. Therefore, we here show some very qualitative discussions on the color superconducting.
To discuss the color superconducting phase with the canonical ensemble method, we must introduce the external field (J ∆ ) and take the J ∆ → 0 limit even at finite θ. Of course, at least for CFL, we can consider the gauge-invariant (color singlet) order-parameter, but it is not possible for the 2SC because unbroken global symmetries in the 2SC phase are the same with those in the hadron phase; for example see Refs. [58,59]. Therefore, we here consider the diquark condensate by introducing the gauge symmetry breaking external field. With the external field, we have the following relation where ∆ means the particular diquark operator; where ∆ C means the expectation value of ∆ operator calculated from the θ region. We here assume that all flavors are degenerated. Since we introduce J ∆ and consider the J ∆ → 0 limit, we can assume where n is infinitesimal value and approaches to 0 in the J ∆ → 0 limit. It should be noted that we may have finite ∆ C (N c n) even if we take J ∆ → 0 because we here consider the sufficiently large but finite size system. Therefore, more strictly speaking, we have ∆ C (N c n) → V where V is the infinitesimal value which approaches to 0 with the V → ∞ limit. In below, we concentrate on the infinitesimal value induced by nonzero J ∆ and thus we simply write V = 0 because following discussions are almost unchanged if we remain V : At least for discussions in this section, the thermodynamic limit is not the matter and thus one can consider following discussions are done with the thermodynamic limit.
As the same discussion in Sec. III and IV, we can expect the exponential suppression factor depends on M B = N c M which can be factored out from the canonical partition function. Therefore, we may write the conden-sate as where f n are coefficients depend on T , n, M and J ∆ at moderate µ R , and the suppression factor is factored out from them; this means that we factor out the suppression factor from Z C (N C n) in the second line of Eq. (32). It should be noted that we assume that exp(−nN c M/T ) factor appears for each canonical sector as in the last line of Eq. (32) because they are nN c -quarks contributions.
In the present work, we are working with T 1 and thus exponential factors can lead the strong suppression at small µ R .
Our interest is that above expression can lead the following situation; ∆ (µ R ) = 0 (low µ R ) nonzero (moderate and high µ R ) .
In other wards, the present question is "Is it possible to manifest Eq. (33) in the canonical ensemble method?".
The answer of the question can be easily understood from the final expression in Eq. (32); if exp[N c (µ R −M )/T ] becomes −1 -order, above situation (33) is realized. Therefore, if the second term is the first nontrivial contribution to ∆ , the necessary energy-scale Λ ∆ for the phase boundary of the color superconducting phase is about µ R ∼ M . When µ R becomes larger value than Λ ∆ , higher-order terms will join the game and thus this discussion will be modified. Detailed value of Λ ∆ is strongly affected by the actual dependence of f n on n; since the higher oscillating modes are expected to be 1/n suppressed with increasing n and thus the second term in Eq. (32) is expected to be the first nontrivial contribution for the nonzero diquark condensate at moderate µ R .
Therefore, the answer for our question "Is it possible to manifest Eq. (33) in the canonical ensemble method?" is yes at least with the present simplified setup. This is the one of possible scenario for the color superconducting from the viewpoint of the canonical sectors.
To discuss one another possible scenario for the situation (33), we consider following setup. To simplify our discussion, we here assume where f is a constant which is less than 1 and denotes the infinitesimal value; f n mimics the 1/n suppression of higher-order oscillating modes. With this assumption, Eq. (32) becomes the sum of infinite geometric series as with r = f e Nc(µR−M )/T .
From the equation, we have the critical chemical potential (µ ∆ ) for the color superconducting as from (1 − r) = . Since f is smaller than 1 in the present setting, µ ∆ is shifted to larger value comparing with M .
To more correct discussion, we need the exact expression of the sum; the coefficient f should have µ R -dependence via Z GC (µ R ) and should remove the divergence. With increasing µ R , the denominator of Eq. (35) can become 1 -order around µ R ∼ M ; its detailed value is affected by the actual value of f . Therefore, we can have the same result obtained in the previous paragraph in the present setup. These means that the situation for the chiral condensate and the diquark condensate are different even if we base on the same canonical ensemble. Above discussions are just theoretical expectation and thus it is not exact because we do not known actual behavior of the second term in Eq. (32) with J ∆ → 0 limit at low T . It should be noted that each n should depends on 1/N c and thus it becomes zero in the large N c limit; it is consistent with the large N c QCD picture. These results can be checked when J ∆ is introduced to the NJLtype model in the canonical ensemble method with the multi-precision calculation used in Refs. [35,36].

VII. SUMMARY
In this paper, we have discussed the confinementdeconfinement nature of QCD at moderate real chemical potential (µ R ) with low temperature (T ) from the viewpoint of the canonical ensembles. The canonical partition function (Z C ) with fixed quark number n is constructed by using the grand-canonical partition function (Z GC ) at finite θ via the Fourier transformation and the fugacity expansion where θ := µ I /T with the imaginary chemical potential (µ I ).
We have assumed that Z GC (θ) is consist of the constant part and the 2π/N c -periodic oscillating mode which is proportional to cos(N c θ); it can be justified at sufficiently low T . With the matching with the energy scale of the quarkyonic phase transition known in the large N c QCD, we can discuss the confinement-deconfinement nature from the change of relevant canonical sectors with varying µ R . It should be noted that this behavior can be obtained in the QCD effective model with some ansatz. There is the balanced region where two canonical sectors are equally contributes to Z GC (µ R ), and it is very tiny in the large N c because the region is 1/N c suppressed.
In the realistic QCD case, it is difficult to make the clear discussion about the quarkyonic phase, but we can still expect that the change of relevant canonical sectors with varying µ R . In this case, the balanced region between nearest-neighbour canonical-sectors is not strongly suppressed and thus there are continues change of the system. This behavior may be consistent with the recent expectation at moderate µ R based on the quakyonic phase and the soft surface delocalization scenario. Unfortunately, the numerical calculation to construct the canonical sectors is quite hard; we need the multi-precision calculation in the Fourier transformation [26,35,36] and thus actual numerical calculation is our future work.
In addition to the quakyonic phase, we have presented some qualitative discussions about the color superconducting and the chiral symmetry restoration from the canonical sectors. This paper guarantees that the canonical ensemble method can work at moderate µ R even with low T and then the quakyonic and color superconducting phases and also the chiral symmetry restoration can be described in the method. Since the exploration of the QCD phase diagram from the lattice QCD simulation with the canonical ensemble method at moderate µ R and low T must be the long journey, we hope present results become the marker for the right direction.