Imaginary Chemical Potential , NJL-Type Model and Confinement – Deconfinement Transition

In this review, we present of an overview of several interesting properties of QCD at finite imaginary chemical potential and those applications to exploring the QCD phase diagram. The most important properties of QCD at a finite imaginary chemical potential are the Roberge–Weiss periodicity and the transition. We summarize how these properties play a crucial role in understanding QCD properties at finite temperature and density. This review covers several topics in the investigation of the QCD phase diagram based on the imaginary chemical potential.


Introduction
Understanding nonperturbative properties of quantum chromodynamics (QCD) at finite temperature (T) and chemical potential (µ) is one of the important and interesting subjects in the elementary particle, nuclear and hadron physics.The confinement-deconfinement transition and the chiral phase transition are famous examples of the properties.To investigate the non-perturbative features of QCD, the lattice QCD simulation is the powerful and gauge invariant approach, but it has the well-known sign problem at finite real chemical potential, µ = (µ R , 0), and thus we cannot perform it exactly in the region.Several methods are proposed so far to circumvent the sign problem, see Ref. [1] as an example, but the applicable regions are still limited in µ R /T < 1.It should be noted that there are some attempts to overcome the limitation: famous methods are the complex Langevin method [2,3], the Lefschetz thimble method [4-6], the path optimization method [7,8] and so on.However, these methods still have serious problems.
Because of the sign problem, low energy effective models of QCD are widely used to explore the QCD phase diagram.The Nambu-Jona-Lasinio (NJL) model is one of the famous low energy effective models of QCD; it can describe the spontaneous chiral symmetry breaking and also its restoration.However, this model cannot describe another important nonperturbative property of QCD, i.e., the confinement phenomena.In Ref. [9], however, the author introduces the Polyakov loop which becomes the order parameter of the confinement-deconfinement transition in the pure Yang-Mills theory to the NJL model and then the door for discussions of the confinement-deconfinement transition in the QCD effective model has been opened; this model, which can describe the chiral phase transition and the approximated confinement-deconfinement transition, is called the Polyakov-loop extended NJL (PNJL) model.It should be noted that there are some other approaches which can treat these properties such as the Dyson-Schwinger equation, the functional renormalization group and the covariant spectator theory; see Refs.[10][11][12][13].
The Polyakov loop which is related with the gauge invariant holonomy is the order parameter of the spontaneous Z 3 (center) symmetry breaking in the pure Yang-Mills theory.The definition of the Polyakov loop is where N c is the number of color, β does the inverse temperature (1/T) and P represents the path ordering.The Polyakov loop can be expressed in the pure Yang-Mills theory as Φ ∼ e −βF q , where F q means the one-quark excitation free energy.Therefore, Φ = 0 means the confined phase because we need the infinite F q to pick up the single quark from the thermal system.In comparison, we can pick up the single quark from the system when Φ is nonzero; it is nothing but the deconfined phase.
Of course, the Z 3 symmetry is explicitly broken when the dynamical quarks are taken into account in the system and then the relationship between the Polyakov loop and the free energy becomes unclear.However, we are still expecting that the Polyakov loop can represent some parts of the confinement phenomena in the dynamical quark system.Actually, the PNJL model with some modifications can well reproduce states-of-arts lattice QCD data at finite T; for example, see Ref. [14] as an example and also see Refs.[15][16][17][18][19] for some lattice data.Although effective models of QCD can provide us several important and interesting pieces of knowledge, such effective models have large model ambiguities.Thus, quantitative predictions are impossible at present.Particularly, ambiguities become serious at moderate and high density regions because we cannot access the lattice QCD data in the region because of the sign problem.This indicates that we need some more approaches to remove model ambiguities to quantitatively discuss nonperturbative natures of QCD at finite density when we employ the effective models of QCD.The imaginary chemical potential, µ = (0, µ I ), which may cause readers to feel the strange sound, but this region is the promising laboratory to remove the ambiguities of the effective models of QCD.This region has the following desirable properties: 1.
Possible analytic continuation process to the real chemical potential, 3.
Relationship with the real chemical potential region via the canonical ensemble.
These properties open the door to remove the model ambiguities, particularly ambiguities induced from the chemical potential effects.The importance of the imaginary chemical potential region can be understood from the canonical ensemble method [20][21][22][23][24] which is related with list 3; the grand-canonical partition function with µ = (0, µ R ) can be obtained from the grand-canonical partition function with µ = (0, µ I ) via the Fourier transformation and the fugacity expansion [25].This fact means that the imaginary chemical potential region has almost all of information of the real chemical potential region.Details are explained later.
In addition to the sign problem, the imaginary chemical potential may provide us with important and interesting knowledge about the confinement phenomena: in ordinary understanding of the confinement-deconfinement transition at a finite temperature in QCD, there are no "phase transitions", i.e., the crossover; there are no singularities in local order-parameters and also the thermodynamic quantities.However, it has been recently discussed in Ref. [26] that the confinement and deconfinement states at zero temperature can be clarified via the topological order [27] and then it is not necessary that the local order-parameter exists and several observable show singular behaviors.The analogy of the topological order has been applied to the thermal QCD by employing the imaginary chemical potential and then it is expected that confinement-deconfinement transition may be determined from the topological viewpoint [28][29][30].This suggests that the imaginary chemical potential can give us important information of the confinement-deconfinement transition.
This review is organized as follows.In the next section, we explain several important properties of QCD at finite imaginary chemical potential such as the Roberge-Weiss periodicity and the transition.
Section 3 discuss what the imaginary chemical potential is.Properties of the NJL-type models are discussed in Section 4. Applications of the imaginary chemical potential to the investigation of the QCD phase diagram are shown in Sections 5 and 6.Section 7 is devoted to Conclusions.

Roberge-Weiss Periodicity and Transition
At finite imaginary chemical potential, there is a special periodicity that is so called Roberge-Weiss (RW) periodicity and a special first-order transition, which is called the RW transition [25].This periodicity was a model independently proven from the structure of the QCD partition function.The special first-order transition was predicted from the strong coupling QCD and the perturbative one-loop effective potential with a background gauge field.It should be noted that these properties are already checked by using the lattice QCD simulations [31][32][33][34][35][36][37][38][39][40].
For reader's convenience, we here show the brief glossary for important features of QCD at a finite imaginary chemical potential, which are the main issues in this review:

•
Roberge-Weiss (RW) periodicity: Special 2πk/N C periodicity of the grand-canonical partition function (Z) along the θ-axis where Symmetry which characterizes the RW transition line.This symmetry is associated from the time reversal or the charge conjugation and Z N c transformations via the semidirect product (It is first discussed by using the combination of the charge conjugation and Z N c symmetries in Ref. [41].).In other words, the system symmetry at θ = (2k − 1)π/N c is enhanced.The modified Polyakov-loop then becomes the order-parameter of the spontaneous breaking of this symmetry.See Section 2.4 for details.

• Roberge-Weiss (RW) endpoint:
Endpoint of the first-order RW transition line.There are possibilities that the endpoint becomes the second-order (trivial scenario) or the first-order (nontrivial scenario) near the physical quark mass.The RW endpoint temperature is denoted by T RW .See Section 2.3 for details.
If readers know the above properties well, the following section can be skipped.

Roberge-Weiss Periodicity
The QCD grand-canonical partition function is where q is the quark field with N f flavor and N c color, γ means the gamma matrices, D represents the covariant derivative which includes the gauge field A µ , g means the gauge coupling and F µν the field strength where some indices are omitted.The gauge fixing and ghost terms are neglected here because they do not affect the following discussions.When we redefine the quark field as q(τ, x) → e iτθT q(τ, x), q(β, x) = −e iθ q(0, x), the chemical potential disappears from the expression of the grand-canonical partition function, but it appears in the temporal boundary condition of the quark field.After performing the (temporal) twisted gauge transformation (U): with where z k := exp(i 2πk N c ), the quark field feels the transformation as This expression of the boundary condition brings us to the fact that: This periodicity is the RW periodicity and the above explanation is the model independent proof.It should be noted that the gluonic part is invariant under the twisted gauge transformation; it is the origin of the Z N c symmetry in the pure gauge limit.In other words, the remnant of the Z N c symmetry in the pure Yang-Mills theory appeared as the RW periodicity in the full QCD.

Roberge-Weiss Transition
Interestingly, behaviors of the RW periodicity are quite different in the low and high T regions.Then, we can find the special first-order transition along T-axis at θ = (2k − 1)π/N c with sufficiently large T which is the so called Roberge-Weiss transition.Unfortunately, the RW transition is difficult to discuss the model independently and thus we here employ the perturbative calculation and the strong coupling limit of QCD as in Ref. [25].It should be noted that this transition is confirmed from the lattice QCD simulations [31][32][33][34][35] and also computations of effective models of QCD [42].
From the following two subsections, we can understand that the Z N c images and the baryonic fugacity play a crucial role for the RW periodicity at high and low T regions, respectively.

RW Periodicity in the Confined Phase
To briefly discuss the confined phase which appears at low T, we here consider the strong coupling limit of QCD with the mean-field approximation [43,44] because it leads us to the confined thermal system.The effective potential of QCD in the strong coupling limit becomes where we only pick up µ-dependent terms because µ-independent terms cannot affect the RW periodicity.Here, we use the analytic continuation from the real to imaginary chemical potentials.
The RW periodicity is induced by the factor cos(N c θ) in Equation (8); this factor comes from the baryonic fugacity; We never obtain the RW transition in this case since there is no origin of the singularity at θ = (2k − 1)π/N c .This fact suggests that the absence of the RW transition is a necessary condition for the realization of the confined phase because we need sufficient strength of the single or some color non-singlet combinations of quark excitations to make the singularity.These results can also be obtained from hadronic models: for example, if we use the chiral perturbation theory with the relativistic Virial expansion [45] and that with the finite energy sum rule [46], we cannot obtain the RW transition because the RW periodicity is induced by the baryon fugacity and thus there is no origin of the singularity.Since we cannot have singularities at θ = (2k − 1)π/3, θ-odd quantities such as the quark number density cannot have nonzero value at θ = (2k − 1)π/3 and then the phase of the Polyakov-loop which plays a crucial role to clarify the Z 3 images at high T are always smooth along the θ-axis.In other words, there is only one Z 3 image at low T.

RW Periodicity in the Deconfined Phase
The deconfined phase can be discussed by using the perturbative calculation and thus we here consider the perturbative one-loop effective potential with the background gauge field [47,48].The RW periodicity in the perturbative one-loop effective potential with the background gauge field is induced by the combination of gA 4 /T + θ; it contains the single quark excitation and it is related with the quark fugacity ξ q = e ±i(gA 4 /T+θ) . ( Below, we set N c = 3 because we are interested in QCD.
The perturbative one-loop effective potential was obtained in Refs.[47,48] and the actual form is with where m is the quark mass, q i is defined as A 4 = 2πTφ with φ = diag(q 1 , q 2 , q 3 ), q ij = q i − q j and K 2 (x) means the modified Bessel function of the second kind.The first and second lines of Equation ( 12) represent the quark and gluon contributions, respectively.The summation over n in the gluonic effective potential can be performed as where B i (x) is the i-th order Bernoulli polynomial with 0 ≤ x ≤ 1.This effective potential can describe the deconfined phase and can be easily extended to arbitral dimensional systems [49].Here, we use the Polyakov-gauge, ∂ 4 A 4 = 0, and remaining gauge freedoms of the space components is used to make the background gauge field diagonal.
In the present case, we find the cusp of the effective potential and the gap of the quark number density at θ = (2k − 1)π/3 and then we can find the first-order RW transition; for example, see Figure 1 which is the PNJL model result, but the PNJL model provides us the same behavior of the perturbative effective potential at high T. (For example, see Refs.[42,50] for results of the PNJL model) and Ref. [25] for results of the perturbative one-loop effective potential.It is well known that θ-even quantities such as that the effective potential, the chiral condensate and the entropy density should have the cusp and θ-odd quantities such as the quark number density should have the gap along the θ-axis.This fact can be understood from the co-existence theorem [51].It should be noted that the singularity is induced by the existence of the Z 3 images: It is well known that there are Z N c images if there is the RW transition along the θ-axis.These images are corresponding to a minimum of the thermodynamic potential.One of them should be the global minima in a certain range of θ, but different minima which is the local minima of the thermodynamic potential in the range can be changed into the global minima in another range of θ.The images can be characterized by the phase of the Polyakov loop: Actually, are all possible Z N c images in the case of N c = 3 at sufficiently high T.

Roberge-Weiss Endpoint
The RW endpoint is the endpoint of the first-order RW transition because the RW transition appears at high T, but it cannot appear at low T.This means that there should be the endpoint at certain T, which is usually denoted by T RW .Since the singularity of the RW transition is the first-order (gap), there is the possibility that the endpoint is the first-order or the second-order.In the heavy quark-mass (m) region, the spontaneous Z 3 symmetry breaking still survives and thus the RW endpoint should be the first-order and then the endpoint becomes the triple-point where three first-order lines meet.When the quark mass becomes lighter and lighter, the endpoint starts to show the second-order behavior from certain m c .However, near the physical quark mass regime, the order of the endpoint may depend on the number of flavor.Some lattice QCD simulations indicate that the order of the RW endpoint is changed into the first-order again and then the RW endpoint becomes the triple-point in the two-flavor and three-flavor systems [36,52,53]; there should be another m c in the system.However, some lattice QCD simulations indicate the second-order RW endpoint even below the physical quark mass in the 2 + 1 flavor system [54][55][56].

Shift Symmetry Breaking
To discuss properties of QCD at finite imaginary chemical potential, the modified Polyakov-loop (Ψ) is very useful.The definition of Ψ and it conjugate are On the RW transition line at θ = (2k − 1)π/3, the quark number density and the imaginary part of the modified Polyakov-loop have the nonzero value above T RW [42].This transition can be understood from the spontaneous shift symmetry ((Z 2 ) shift ) breaking.The (Z 2 ) shift symmetry is the invariance under the Z 2 transformation associated from the time reversal (T ) or the charge conjugation (C) and Z 3 transformations via the semidirect product [41,[57][58][59], i.e., Z 2 Z 3 .
Actual behavior of the Polyakov loop and the modified Polyakov-loop under the (Z 2 ) shift transformation becomes at θ = π/3, but the imaginary part of the modified Polyakov-loop is transformed as where A − −−− → B × C means that A is not transformed to C by using B transformation.Therefore, we can use ImΨ as the order parameter to detect the spontaneous (Z 2 ) shift symmetry breaking at θ = π/3.

Interplay of Imaginary Chemical Potential
It is well known that the imaginary chemical potential region is related with the real chemical potential region via the Fourier transformation and the fugacity expansion [25]; where Z C is the canonical partition function, Z GC means the grand-canonical partition function and Q is the real quark number.This means that the imaginary chemical potential region has the information of the real chemical potential region because the most basic quantity, Z GC (T, µ R ), in the thermal system can be exactly constructed from Z (T, µ I ).Therefore, if we can construct the reliable effective model at finite µ I , it means that such model is also reliable at finite µ R .The error-bar of lattice QCD data affects the error of the effective model at large µ R , but it is controllable in principle.These procedures for the model building of QCD is the so-called imaginary chemical potential matching approach [60].The actual flow chart of the imaginary chemical potential matching approach is as follows: 1.
Prepare lattice QCD data for several observables at finite θ.

2.
Prepare a suitable effective model which reproduces the RW periodicity and the transition.

4.
Calculate observables by using the model and compare them.
Until we can have a good parameter set that well reproduces the lattice QCD data, we repeat steps 4 and 5.
From the above explanations, the importance of the imaginary chemical potential can be understood.However, the meaning of the imaginary chemical potential is still unclear because it leads us to the pure imaginary quark number density: at finite µ I , we can define the quark number density as Thus, it becomes the pure imaginary and then it seems unphysical.Of course, there should be discussions that this quantity can be interpreted as the quark number density or not because we can interplay it as different valuables as explained below.

Analytic Continuation
Simplest interpretation of the imaginary chemical potential is the analytic continued value of the real chemical potential.Then, the negative and the positive region of µ 2 are corresponding to the imaginary and real chemical potentials, respectively.Based on this interpretation, several investigations of the QCD phase diagram have been done via the analytic continuation method.Details of the method is explained in Section 5.1.

Boundary Condition of Fermion for the Temporal Direction
There are some more different interpretations of the imaginary chemical potential.The famous one is that the imaginary chemical potential can be converted to the fermion boundary condition for the temporal direction.This can be easily seen from the Matsubara frequency with the arbitral boundary condition as where ω n is the usual fermion Matsubara frequency, and φ is the boundary angle and φ = φ − 1/2.Therefore, factor 2πTφ can be redefined as µ I .In this case, ∂V /∂µ does not have any meaning and then we can be free from the interpretation problem of the pure imaginary quark number density.This fermion boundary condition can play an important role in the adjoint QCD and also in gauge-Higgs unification models because it causes the Hosotani mechanism.In addition, we can construct the Polyakov-loop like quantity by considering the boundary-angle dependent chiral-condensate.This quantity is the so-called dual quark condensate [61] and its definition is where n means the winding number.The Polyakov loop is the winding number 1 and thus Σ (1) behaves similarly to it.Since the dual quark condensate is constructed from the chiral condensate, this quantity is a good quantity to investigate the correlation between the chiral and deconfinement transitions.
The interesting consequence from the dual quark condensate is that it increases with increasing T also in the NJL model [71][72][73].The NJL model is the Φ = 1 limit of the PNJL model and thus this behavior is strange because there is no spontaneous Z 3 symmetry breaking contribution.This strange behavior in the NJL model was clarified in Ref. [74].The behavior of the dual quark condensate in NJL model is coming from the chiral transition because the quark mass that is affected by the chiral symmetry breaking directly affects the explicit Z 3 symmetry breaking.The same effects coming from the chiral transition should be there in the Polyakov loop and thus different determination of the deconfinement transition without using the Polyakov loop is interesting.

Aharonov-Bohm Phase
The other interpretation is the Aharonov-Bohm phase [75].At a finite temperature, the temporal direction is compactified in the imaginary time formalism and then we can insert U(1) flux, θ, to the (fictitious) hole.Because of the Aharonov-Bohm effect, the vector potential a 4 (θ) = Tθ is induced and the appearance form of a 4 in the action is exactly the same as the imaginary chemical potential; for example, see Ref. [76].In this case, ∂V /∂θ becomes the surface current density.In the imaginary time formalism, the time direction becomes the imaginary time direction to introduce the temperature to the system and thus such surface current density is difficult to understand.It may be related with the problem that the Z 3 images can be determined in the Minkowski space or not; it remains an open question.
If we obey the above interpretation, we can consider following this interesting scenario.The interpretation from the Aharonov-Bohm phase may have one advantage comparing with some other interpretations because we can use the discussion of the topological order [27].Applications of the topological order to zero temperature QCD was discussed in Ref. [26].In Ref. [26], the authors consider the torus T 3 at zero temperature.There are three important adiabatic operations: 1.
U(1) flux insertion to holes of spatial closed loops.

3.
Moving of the quark along loops.
These operations are represented by U a , σ i and τ a i , where a represents the space indices.Commutation relations of σ i and τ a i are described by the Braid group, and the commutation relations of those with U a are characterized by the Aharonov-Bohm effect.If the fundamental degree of freedom in the system is quark, the commutation relation becomes non-commutable because of the fractional charge, but it is commutable if hadrons are fundamental degree of freedom in the system.Therefore, if there is only one vacuum in the former case, it is inconsistent with the non-commutability of the operations and thus vacuum degeneracy should be there.The actual degeneracy in the deconfined phase is the 3 3 for the case of T 3 tours system.On the other hand, the confined phase does not have such degeneracy.
In the discussion of Ref. [26], they consider T 3 tours and obtain the nontrivial vacuum degeneracy for the deconfined state.On the other hand, we can consider U(1) flux insertion to the hole of the compactified time loop: we are considering the finite T system, but the topological order cannot be well defined in the system because the finite T state is constructed by a superposition of zero T states with a Boltzmann factor.Therefore, we cannot operate above the operations adiabatically.However, the RW periodicity shows the significant difference between them as already explained in Sections 2.2.1 and 2.2.2.It is induced by the nontrivial appearance of the RW periodicity in the deconfined phase.This may mean that we can clarify the confinement and deconfinement phases from the non-trivial degeneracy of the effective potential and it seems similar with the vacuum degeneracy in the zero temperature system.Thus, we may distinguish two phases from the structure of the periodicity.We may then state that the RW endpoint temperature is the critical temperature of the deconfinement transition; of course, it is still a conjecture.Actually, this definition is perfectly consistent with the definition of the critical temperature of the deconfinement transition by using the susceptibility of the Polyakov loop in the infinite quark mass limit where Polyakov loop is the exact order-parameter of the deconfinement transition.
In the discussion of the topological order, the phase transition can not be described by the order parameter, which is associated with the spontaneous symmetry breaking.At θ = (2k − 1)π/3, there is the singularity and then the quark number density can has nonzero value above T RW .This suggests that the imaginary part of the quark number density at θ = (2k − 1)π/3 can be used as the order parameter like a quantity of the deconfinement transition in the above scenario.In Ref. [29], the quark number holonomy has been proposed based on the above discussions.The functional form of the quark number holonomy is defined as where ñq is a dimensionless quark number density such as n q /T 3 .The quark number holonomy can count gapped points of the quark number density in the −π ≤ θ ≤ π region.Therefore, it is characterized by the nontrivial free-energy degeneracy because the quark number density should have the gap at θ = (2k − 1)π/3 when the free-energy is non-trivially degenerated.

NJL-Type Model at Finite Imaginary Chemical Potential
To discuss the finite density QCD, we cannot perform lattice QCD simulation at present because of the sign problem.Thus, we should employ the effective model.Below, we explain the NJL-type model as the typical QCD effective model that can be used at finite density.The system is set to the two flavor and three color case.

Nambu-Jona-Lasinio Model
Our starting effective model is the Nambu-Jona-Lasinio (NJL) model and the Lagrangian density of the two-flavor and three-color NJL model in the Euclidean space is where m 0 means the current quark mass and the constant G is the coupling constant of the scalar and pseudo-scalar type four-Fermi interactions.This model is the simplest low-energy effective model of QCD, but this model only has the trivial 2π periodicity for θ-direction.This fact can be easily seen from the effective potential with the mean-field approximation as where The condensate σ is defined as σ ≡ qq .Because of the exponential factor of the effective potential, we can find the trivial 2π periodicity when we replace µ with iµ I .
On the other hand, the Polyakov-loop extended NJL (PNJL) model [9] can overcome this problem as explained in the next subsection.The foundation of the NJL model and also the PNJL model are shown in Ref. [77,78].In the NJL model and the PNJL model, we can consider some more interactions such as a vector-type interaction.Details of effects of the vector-type interaction which are important for the finite chemical potential have been investigated in the NJL model [79,80] and the PNJL model [81].Determination of this interaction is very difficult in the usual model construction because the NJL-model parameters are determined at zero chemical potential, but this interaction becomes the density-density interaction in the mean-field approximation.One method to determine it is using the small µ R information for the quark number density of lattice QCD data [82].In comparison, the strength of the vector-type interaction can be determined directly via the quark number density in the imaginary chemical potential matching approach; for example, see Ref. [42].

Polyakov-Loop Extended Nambu-Jona-Lasinio Model
The PNJL model can approximately describe the deconfinement and chiral phase transitions at the same time by introducing non-perturbative effects through the NJL model and Polyakov-loop effective potential [9].The Lagrangian density of the two-flavor and three-color PNJL model in the Euclidean space is where the covariant derivative is D ν = ∂ ν − igA ν δ ν4 and V PNJL,g expresses the gluonic contribution.
The actual form of the effective potential with the mean-field approximation is where with the Fermi-Dirac distribution functions; where Φ and Φ are the Polyakov loop and its conjugate.In the PNJL model, the Polyakov loop is In the correct definition of the Polyakov-loop, we should take the expectation value of tr c e iA 4 /T itself, but this definition is difficult in the PNL model and thus we usually use a gauge variant expectation value to define it.From the Jensen inequality, the model Polyakov-loop provides the upper-bound of the correct Polyakov-loop [83] as tr e igA 4 /T ≤ tr e ig A 4 /T .Therefore, we can use it if we suitably set model parameters to reproduce a scale of the deconfinement transition in the pure gauge limit.Here, we use the logarithmic Polyakov-loop effective potential [84] as the effective model for the gluonic contributions.The functional form is with where parameters, (a 0 , a 1 , a 2 , b 3 ), are shown in Ref. [84].The remaining parameter T 0 is usually fixed as 270 MeV, which is the critical temperature in the pure gauge limit.In the case of full QCD, T 0 is determined to reproduce the pseudo-critical temperature of the deconfinement transition at zero µ; there should be the back-reaction of quark contributions to the Polyakov-loop effective potential [85,86].
Here, we use T 0 = 270 MeV because we are interested in the qualitative behavior.We here use the same gauge fixing and path integral formulation of the perturbative one-loop effective potential.In addition, parameters of the NJL part used in this study is employed from Ref. [81].
It should be noted that we can see the importance of the modified Polyakov-loop in the PNJL model because the Fermi-Dirac distribution functions of the PNJL model can be expressed as and thus the effective potential of the PNJL model clearly has the RW periodicity because the gluonic part is invariant under the Z 3 transformation and the modified Polyakov-loop and the factor exp(∓3iθ) have the RW periodicity by definition.For reader's convenience, we here show the θ-dependence of the chiral condensate σ/σ 0 and the normalized quark number density n q /T 3 in Figure 1.Since the chiral condensate (the quark number density) can have the θ-even (θ-odd) quantity, it can have the cusp (gap) at θ = (2k − 1)π/3 in the case with sufficiently large T. The chiral condensate is normalized by using σ at T = µ = 0. We can clearly see the RW transition at T = 300 MeV and then the chiral condensate has the cusp and the quark number density has the gap at θ = (2k − 1)π/3, but not at low T. Therefore, the PNJL model can reproduce the RW periodicity, RW transition and RW endpoint, automatically and systematically.
Other effective models for the gluonic contributions are the Meisinger-Miller-Ogilvie model [87], the Polynomial-type Polyakov-loop effective model [88][89][90][91], the matrix model for deconfinement transition [92] and the effective potential from the Landau-gauge gluon and ghost propagators [93].When we suitably set parameters of the models, those effective models can reproduce lattice thermodynamics at zero chemical potential in the same level in principle.It is, however, changed when we consider the heavy quark mass region and we can check and remove the model ambiguity of the gluonic contribution by considering the Columbia plot [94,95].Particularly, the heavy quark mass region, where critical quark mass (m c ) appears, can be used as the laboratory for this purpose because the actual value of m c strongly depends on the model details [95]; see Section 2.3 for the explanation of m c .

Application of Imaginary Chemical Potential to Explore the QCD Phase Diagram
There are several ways to utilize the imaginary chemical potential to investigate real chemical potential region.Here, we explain three different famous methods, the analytic continuation method, the canonical ensemble method and the Lee-Yang zero analysis as an example.In the explanation of the canonical ensemble method, we also discuss the Polyakov-loop paradox which is the serious problem because it makes the foundation of the method unclear.

Analytic Continuation Method
In the analytic continuation method, we first generate lattice configurations at each θ in the range 0 ≤ θ < π/3.Actually, there is no sign problem and thus we can exactly perform the lattice QCD simulation.After generating the configurations, we calculate some observables which we wish to know at finite µ R .The obtained results are fitted by using analytic functions such as the Polynomial function (F) and then we can analytically continue the quantity to the positive µ 2 region from the negative µ 2 region.Of course, some other functional forms can be used such as the trigonometric functions and the Pade approximation.In some observables, simple oscillation functions based on the Fourier decomposition can work well to fit the lattice QCD data; see Ref. [31,34,96] as an example.
It should be noted that thermodynamic quantities and order parameters are periodic in terms of µ I /T.Unfortunately, we have the first-order RW transition at θ = π/3 where the singularity appears and thus it leads the convergence radius to the method.Therefore, we cannot go beyond µ R /T = π/3, in principle.In addition, if there are singularities which are induced from the first-order chiral phase transition at finite µ R /T < π/3, we cannot go beyond the point.Therefore, this approach is strongly limited from the singularity issue.To see these problems, the Gross-Neveu model is a good laboratory even if it does not have the RW periodicity [97]; we can understand how higher-order terms contribute to the convergence.

Canonical Ensemble Method
As mentioned above, the imaginary chemical potential has the information of the real chemical potential via the Fourier transformation and the fugacity expansion.Actually, we can construct the canonical partition function Z C with real quark number Q as and then we can have By using these relations, we can indirectly calculate the real chemical potential region by using the lattice QCD simulation.Unfortunately, the numerical error becomes large when we decrease the temperature and thus it is difficult to reach the low T and high µ R region.
Since the QCD grand-canonical partition function has the RW periodicity, Equation ( 32) can be rewritten as where Therefore, Z C (T, Q) becomes Then, the expectation value of an arbitrary RW periodic operator (O) is given by where • • • GC means the expectation value directly evaluated from the grand-canonical ensemble and the term "undefined" means that we encounter 0/0 in the calculation, but it can be defined if we remove the total factor which leads 0 by reducing the fractions to a common denominator because the factor appears in both the denominator and numerator.
The RW periodicity induces the Polyakov-loop paradox in the canonical ensemble method which means that the Polyakov loop becomes exactly zero when we use the canonical ensemble; where we limit the summation of Q by B = 3Q to avoid the divergence because of the RW periodicity.The Polyakov loop is not the RW periodic quantity and thus it does not follow Equation (37); it is the origin of the Polyakov-loop paradox.Then, we always obtain Φ GC = 0 = Φ GC .Of course, some other RW periodic quantities such as the chiral condensate are not affected by the paradox, but the RW un-periodic quantities such as the Polyakov loop feel the paradox.Because of the Polyakov-loop paradox, the foundation of the canonical ensemble method seems to be unclear.However, in Ref. [98], the authors resolved the paradox by considering two different approaches; one is the modified Polyakov-loop representation and the other is the trivial Z 3 -image restriction.

Approach 1: Modified Polyakov-Loop Representation
The modified Polyakov-loop is the RW periodic quantity [42] and thus it can be nonzero in the canonical ensemble method; Since we have which is mathematically true because of properties of the Fourier transformation, the Fourier transformation can be evaluated as Therefore, we should have Ψ = Φ at µ R = 0 by shifting (Q − 1) to Q as when Q and Q run from −∞ to ∞.Therefore, the Polyakov-loop paradox is gone by using the modified Polyakov-loop.This may mean that the modified Polyakov-loop is a more fundamental quantity than the Polyakov loop in the imaginary chemical potential.

Approach 2: Trivial Z 3 -Image Restriction
This paradox can be evaded by using the modified Polyakov-loop or the restriction of the integral range as −π/3 ≤ θ ≤ π/3; It should be noted that this restricted Fourier transformation can exactly reproduce the non-restricted result and this fact is mathematically justified in Ref. [98].Since we do not have the total factor appearing in Equation (36) in this restricted integration because the total factor is induced from the existence of the non-trivial Z 3 -images, these are removed.Therefore, we have From this subsection, we can understand the importance of the RW periodicity and its control because the RW periodicity seriously affects the canonical sectors and then its treatment is crucial, i.e., evading the Polyakov-loop paradox.The imaginary chemical potential is related with the canonical sector and thus it can have lots of information related to the confinement-deconfinement transition because the canonical sectors are expected to have information of the quark excitation modes in the thermal system.

Lee-Yang Zero Analysis
The analytic continuation and canonical ensemble methods use the imaginary chemical potential as the reference system to investigate the real chemical potential region.In comparison, there is a method which more directly uses the imaginary chemical potential to explore the QCD phase diagram; it is a so-called Lee-Yang zero analysis [99].
In the Lee-Yang zero analysis, we introduce the complex chemical potential (or some other complex values such as the complex temperature) and then we search zeros of the partition function on the complex plane.Then, the canonical approach helps us to treat the complex-valued grand-canonical partition function in QCD.Usually, we can perform the simulation in the finite size system, and zeros cannot appear just on the µ R axis because standard phase transitions are smeared by the finite size effect, but these can appear on the complex plane.If we enlarge the system size, zeros appearing near the µ R axis on the complex plane finally appear on the µ R axis in the continuum limit if there are phase transitions.Therefore, we can clarify the existence of the phase transition from the approaching behavior of zeros to the µ R axis.Unfortunately, we do not have clear evidence that there are phase transitions at finite µ R from the Lee-Yang zero analysis yet, but there is clear evidence for the existence of the RW transition at θ = (2k − 1)π/3 from the lattice QCD simulations [100-102].

Similarities Measurement
In this section, we discuss the similarity between the PNJL model and QCD as part of this review.As mentioned above, the imaginary chemical potential region can be used as the laboratory to check which effective models can reproduce the lattice QCD data.In this case, it is important to know how similar the effective model and lattice QCD are.The coefficients of the Fourier decomposition [103,104] are promising quantities for this purpose because such coefficient is related with the higher-order cumulants which are expected to be relevant to understand the confinement-deconfinement transition; for example, see Ref. [105,106].Unfortunately, such Fourier decomposition seems to be numerically difficult at present because we must fit lattice data or the effective model results at the finite imaginary chemical potential, but it needs very accurate numerical data to pick up higher-order components; see Refs.[30,39] as an example.
Instead of using the Fourier decomposition, we here consider the "similarity measurement" which may be another promising quantity for the purpose.When we can convert the quark number density as the probability distribution calculated from the theory and wish to know the similarities with simple oscillation functions which can have the relation with the exciting modes in the thermal system explained later, the similarity measure can play an important role.There are several similarity measures such as the Euclidean distance, Manhattan distance, cosine distance, Jaccard measure, Kullback-Leibler divergence, Jensen-Shannon divergence and so on.Here, we choice the Jensen-Shannon divergence [107] for discrete probability distributions, p and q, defined as because it has the following useful relations; D JS (p||q) = D JS (q||p), 0 ≤ D JS (p||M) ≤ ln 2, where M i = (p i + q i )/2 and N mean the number of data points.The triangle inequality is not manifested which is one of the essential properties of distance and thus it is called divergence.When the probability distributions are perfectly matched with each other, D JS becomes exactly zero.In comparison, D JS becomes ln 2 which is the asymptotic value when there is no overlap between the probability distributions.In the case with the continuous probability distribution, the summation in Equation ( 45) is replaced by the integration.In this section, we wish to pick up the non-trivial quark-gluon dynamics from the oscillating behavior of the quark number density and thus we prepare the probability distribution as a function of θ, p(θ i ) where N p is the normalization functor and n (PNJL) q is the quark number density observed on the lattice.We choice the set of reference probability distributions, q J (θ i ), as N q J is the normalization factor and J = 1, 2, 3, • • • .The functions with J = 1, 2 and 3 can be considered as the single, double and triple quark like contributions.It should be noted that the imaginary chemical potential does not introduce the additional energy scale to the system when we fix T because θ can be transformed into the quark temporal boundary condition and does not create the Fermi surface unlike the real chemical potential.In addition, this way has the following desirable properties; we may be free from the renormalization problem on lattice because the absolute value of the quark number density does not matter.It can also be calculated in the finite size system without caring about the finite size effect unlike the quark number holonomy (21), which is always zero without the careful extrapolation to the continuum [29].Data points appearing near the bottom edge are corresponding to the results in the high temperature phase and the data points appearing near the top vertex are the results in the low temperature phase.In both cases, data are well separated; data for low and high T appear around the top vertex and the center of the bottom edge of the 2-simplex, respectively.It indicates that the PNJL model can well reproduce the oscillating modes in lattice QCD data.Then, we can understand the importance of the RW periodicity and also the transition to understand the confinement-deconfinement phenomena.Therefore, the PNJL model can be considered good element assembly to construct the reliable QCD effective model due to the oscillating modes, which directly affect the reliability of the model at finite real chemical potential, existing in the model.Actually, if theory or effective models do not have such oscillating modes, we cannot discuss higher-order cumulants, in principle.

Conclusions
In this review, we overview QCD properties at finite imaginary chemical potential and its applications to exploring the QCD phase diagram at a finite real chemical potential.In particular, we explain why the Polyakov-loop extended Nambu-Jona-Lasinio (PNJL) model can be considered as the promising low-energy effective model of QCD.
First, we have explained several important and interesting features of QCD at finite imaginary chemical potential such as the Roberge-Weiss (RW) periodicity and the RW transition based on Ref. [25].The RW periodicity explains the model independently.In comparison, we employ the strong coupling limit of QCD and the perturbative one-loop effective potential to explain the RW transition.The origin of the RW periodicity is different in the low and high temperature phases; this difference induces the first-order RW transition where the phase of the Polyakov loop and the quark number density have the gap.To understand the RW transition, the spontaneous shift symmetry braking plays a crucial role and then the modified Polyakov-loop is treated as the order-parameter of the symmetry breaking.
In Section 3, it is shown that there is a possibility that we can approach the confinement-deconfinement nature based on several non-trivial interpretations of the imaginary chemical potential such as the temporal fermion boundary angle and the Aharonov-Bohm phase.Based on the interpretations, we have explained why the dual quark-condensate and the quark number holonomy can be considered as the (quantum) order parameter or the indicator of the confinement-deconfinement transition.In addition, we have shown the strategy of the imaginary chemical potential matching approach which combines the low-energy effective model of QCD and the lattice QCD data at finite imaginary chemical potential.
The possible low-energy QCD effective model which can be used in the imaginary chemical potential region has been discussed in Section 4. We have explained why the simple Nambu-Jona-Lasinio (NJL) model can not reproduce the RW periodicity, while the Polyakov-loop extended Nambu-Jona-Lasinio (PNJL) model can.Because of this fact, the PNJL model is the promising low-energy effective model of QCD.Some theoretical problems and confusing points of the PNJL model are also explained.For example, the reason why the model Polyakov-loop cannot be matched with the correct Polyakov-loop is shown.In addition, the detailed reason why the PNJL model has the RW periodicity is explained from the viewpoint of the modified Polyakov-loop; the Fermi-Dirac distribution function of the PNJL model can be expressed by using the modified Polyakov-loop.
In Section 5, we have summarized some approaches to attack the QCD phase diagram at finite real chemical potential via the imaginary chemical potential region; the analytic continuation method, the canonical ensemble method and Lee-Yang zero analysis are explained.Particularly, we show a detailed discussion on the Polyakov-loop paradox appearing in the canonical ensemble method which is important to ensure the foundation of the method.Then, we can see that the usefulness of the modified Polyakov-loop and also the removal of the nontrivial Z 3 -images from the the canonical ensemble.
Finally, the similarity measure has been discussed to see how the PNJL model and the lattice QCD simulation provide the similar oscillating behavior at finite imaginary chemical potential.Actually, we employ the Jensen-Shannon divergence to measure the similarity.From the results, we can understand that the PNJL model can be used as a promising prime-field to model the QCD properties.
Funding: K.K. is supported by the Grants-in-Aid for Scientific Research from JSPS, No. 18K03618.

Figure 1 .
Figure 1.The θ-dependence of σ/σ 0 and n q /T 3 .The dotted and solid lines represent the result at T = 100 and 300 MeV, respectively.

Figure 2 .
Figure 2. The right and left panels show the D. J JS trajectory on the 2-simplex for the PNJL model and the lattice QCD data [39], respectively.