Thermodynamics of 2+1 Flavor Polyakov-Loop Quark-Meson Model under External Magnetic Field

We study the thermodynamics of QCD system under external magnetic field via the 2+1 flavor Polyakov-loop quark-meson model. To incorporate quantum and thermal fluctuations, the functional renormalization group approach is implemented in our work. Pressure, entropy density, magnetic susceptibility and other thermodynamic quantities are calculated and analyzed to investigate the effect of magnetic field on the QCD system. The calculated results are in reasonable agreement with lattice QCD simulations and perturbation theory. We then give an intuitive picture for the response of QCD system to the magnetic field.


I. INTRODUCTION
Researches on the phase transitions of strong interaction matter (QCD phase transitions) have attracted great attentions in recent years (for example, see Refs. ) since they may shed light on revealing the process and the nature of the early universe evolution. It is known that not only during the early stage of the universe evolution but also in the modern heavy-ion collision experiments, there exists strong magnetic field [37,38], which may influence the property of the strong interaction matter drastically. The effects of the magnetic field have naturally been investigated intensively. As a consequence, it has arisen some interesting questions, for instance, the chiral magnetic effect [39][40][41][42][43][44] and the magnetic catalysis effect [24,[45][46][47][48][49][50].
However, recent lattice QCD simulation results show us that, when considering the interaction more sophisticatedly, the magnetic field does not increase the critical temperature for the chiral symmetry to be restored [25,27,31], which is referred to as the inverse magnetic catalysis (magnetic inhibition) effect. Plenty of theoretical works have then been accomplished to clarify which of the two completely opposite effects is correct, or to what extent each of the effects works [10,11,27,31,, even explore the underlying mechanism of the effects from a fundamental point of view (see, for example, Refs. [45,47,54,57,61,[66][67][68][69][70][71][72][73][74][75][76][77][78][79][80] and recent reviews [59,86,87]). Even though there is still no definitive conclusion, it has been recognized commonly that the external magnetic field produces a spatial dimension reduction for charged particles (see, e.g., Refs. [47,54,76,78,79,83,86]). However one has not yet arrived at a common idea on the effect of the dimension reduction, for instance, Refs. [76,83,86] claim that the dimension reduction strengthens quark condensate and leads to magnetic catalysis at low temperature but * Corresponding author: yxliu@pku.edu.cn Refs. [54,79] say that the low dimensionality may favor the inverse magnetic catalysis. It is then imperative to clarify the effect of the external magnetic field, especially with an intuitive or a pedagogic picture.
It has been well known that the equation of state (EoS) and the entropy density can represent well the properties of a matter. Then the EoS and the entropy density of the strong interaction matter in an external magnetic field have been studied with the lattice QCD simulations [27][28][29], Hadron Resonance Gas (HRG) model [30], quark-meson (QM) model [67,68,88] and Polyakov-loop improved quark-meson (PQM) model [10-13, 59, 89-93]. However, the mechanism of the magnetic field effect on the QCD system has not yet been fully understood by analyzing the thermodynamical properties concretely. For instance, the responses of the quark condensate, pressure and entropy density of QCD system to the magnetic field at different temperature regions and the underlying physics picture remain unclear, which all require more investigations. We then investigate the effect of the magnetic field on the thermodynamic properties of the strong interaction matter and try to give an intuitive view point of the mechanism of the magnetic catalysis effect in some simple models in this paper.
In this work, we take the (2+1) flavor Polyakov-loop quark-meson model [10-13, 59, 89-93] to approximate the QCD system. Compared with previous works [10-12, 59, 67, 68], the effects of Polyakov-loop and strange quarks are taken into considerations. Adding Polyakovloop can imitate the effect of confinement and including strange quark brings more mesons (such as kaon) into the system so that the system is more close to the real one at the high energy scale. Besides using a more realistic model, functional renormalization group (FRG) approach [94][95][96][97] is implemented to investigate the thermodynamics of QCD system at finite temperature and magnetic field. FRG approach can incorporate the fluctuations of quarks and mesons thus goes beyond the meanfield approximation. It has then been shown powerful to solve non-perturbative problems.
The remainder of this paper is organized as follows. In Sec. II, we introduce briefly the (2+1) flavor PQM model. In Sec. III, we depict the framework of the FRG approach. In Sec. IV, we describe our calculation procedures and the obtained results, and provide an intuitive picture for the magnetic effect. In Sec. V, we give our summary and some remarks.

II. 2+1 FLAVOR PQM MODEL
The Polyakov-loop quark-meson model is an effective approximation of QCD system [10-13, 59, 89-93]. The Lagrangian (in Euclidean space) of the (2+1) flavor PQM model under an external magnetic field and at vanishing chemical potential can be simply given as The first line represents quark sector, ψ is the quark field with 3 flavors (u, d, s) and 3 colors. D µ = ∂ µ − iqA em µ is the covariant derivative coupled to quarks with q = diag( 2 3 e, − 1 3 e, − 1 3 e) being a diagonal matrix encoding the charge for each quark. We assume that the external magnetic field B is homogeneous and along the z-direction, the corresponding field A em µ can be written as The Yukawa term gψΣ 5 ψ is also included to represent the interaction between quarks and mesons. The first two terms in the second line of Eq. (1) represent the meson sector. Meson fields can be combined to form matrices, which read where σ a and π a are the scalar and pseudo-scalar meson nonets. T a is the SU(3) generator in flavor space, and The derivative term for mesons is a little complicated which is defined as U (Σ) is the meson potential which contains several parts and reads where U (ρ 1 , ρ 2 ) is the symmetric term under U V (3) × U A (3) transformation. ρ 1 and ρ 2 are two invariants constructed from Σ as h x σ x and h y σ y are the explicit symmetry breaking terms and they are added to reduce the SU V (3) symmetry to SU V (2) which is approximately conserved in QCD. σ x and σ y are related to the meson fields via a rotation ξ is the Kobayashi-Maskawa-'t Hooft determinant which is used to imitate the U A (1) anomaly effect. It breaks the U A (1) symmetry explicitly and reads In order to take color confinement effect into consideration, a background gluon field A 4 and its potential V poly (Φ) are added into the Lagrangian. Φ is called Polyakov-loop and is defined as where P is path order operator. V poly (Φ) is a phenomenological potential obtained by fitting lattice QCD data and it's used to represent the gluon's contribution to the Lagrangian. The specific form of the V poly and the parameters therein will be discussed in Sec. III B.

A. THEORETICAL FRAMEWORK
After the Lagrangian is introduced, several approaches can be taken to calculate the effective action Γ of the system. For example, mean-field approximation is a widely used method. It treats meson fields as background and neglects their fluctuations. However, it has been shown that meson's fluctuations are important for investigating QCD system at low temperature [68]. As a consequence, we take the functional renormalization group (FRG) approach instead in this article. With the FRG approach, one can incorporate the fluctuations of mesons as well as quarks.
FRG approach is usually employed to solve nonperturbative problems (see, e.g., Refs. [94][95][96][97] for details). It is well known that the Γ can be expanded in terms of loops and FRG approach implements a momentum-dependent mass term R k (p) for each propagator in the loops. These modified propagators get highly suppressed at small momentum region (p < k) while are left unchanged at large momentum region (p > k). When k is near the cutoff Λ, all the loop integrals are zero and only bare action contributes to the effective action, i.e., Γ k=Λ = S Bare . When k approaches zero, all the loop integrals are restored and we get the full effective action as Γ Full = Γ k=0 . These ideas can be formulated into a functional differential equation, which reads where Γ 2B k (Γ 2F k ) denotes the functional derivative of Γ k with respect to boson (fermion) field. R B,F k is the momentum-dependent mass term assigned to each particle. Tr represents trace operation in all spaces (momentum, flavor, color, · · · ).
One of the main characters of FRG approach is its differential structure. S Bare serves as the initial condition for Eq. (10). We can study the whole system by solving the differential equations instead of doing loop integrals or solving integral equations.
However, it is extremely hard to solve the functional equations exactly. Some simplifications are necessary to solve the problems practically. We will briefly review the necessary truncations and the derivations for the PQM model as follows. The details can be found in Refs. [59,68,88,[90][91][92][93].
We employ the so called local potential approximation (LPA). In this approximation, all the non-local terms and anomalous dimensions of the fields are omitted. The truncated Γ k now reads where only U k (ρ 1 , ρ 2 ) flows with scale k and other terms are the same as Eq. (1). Substituting Eq. (11) into Eq. (10) and implementing the optimized regulators taken in Refs. [10,59,68,95], one can obtain an analytic flow equation for the U k as For the definitions of α b(f ) , E b(f ) , see Ref. [68]. n b is the ordinary boson distribution function whileñ f is the Polyakov-loop modified fermion distribution function which reads At low temperature, Φ ∼ 0 and Eq. (13) is reduced tõ This distribution implies that the average energy of a quark is nearly one third of an ordinary free fermion. In other words, quarks are statistically confined at low temperature. At high temperature, Φ ∼ 1 andñ f is reduced to the ordinary fermion distribution function.
There are usually two algorithms to solve Eq. (12). One is the so called grid method. It discretizes the arguments of U k (i.e., the ρ 1 , ρ 2 and Φ mentioned above) on a multi-dimensional grid, then Eq. (12) can be transformed into a set of coupled ordinary differential equations. However, grid method requires usually much numerical effort, especially when the number of the arguments increases. So Taylor expansion method can be employed instead. Taylor expansion method expands the U k (ρ 1 , ρ 2 ) around the minimums ρ 1 (k), ρ 2 (k) and the series is truncated to a certain order. The expansion series reads Setting n = 5 has been shown to be sufficient for numerical convergence [93]. A set of coupled ordinary differential equations for a i,j and ρ i (k) can then be obtained by substituting Eq. (15) into Eq. (12). The flow equations read simply (specific form can be seen in Refs. [59,68]) and there are many numerical methods to solve these equations.
In short, effective action Γ k can be parameterized into a set of running coupling constants and the nonperturbative β function for each coupling constant can be obtained by using FRG approach. Full effective action can be obtained as Eq. (16) is solved.

B. NUMERICAL IMPLEMENTATION
In this article we take the Taylor expansion method to carry out the numerical calculations. Solving Eq. (16) requires initial conditions which are the parameters in the Lagrangian, i.e., Eq. (1). We make use of the parameters given in Ref. [98] which are summarized in Table I. g hx/Λ 3 hy/Λ 3 ca/Λ a10/Λ 2 a20 a01 Λ/GeV 6.5 0.121 3 0.336 3 4.808 0.56 2 26 50 1 The first three parameters in Table I are fixed by the Goldberger-Treiman relation and the PCAC theorem [99]. c a represents the strength of U A (1) anomaly and it depends generally on the external parameters such as temperature or magnetic field. In this article we set c a to a constant to simplify the discussions. Λ is the ultraviolet cutoff of PQM model and it's fixed by the result that gluon develops a mass gap below 1 GeV and then decouples with low energy physics [100]. Thus in practical calculations, we tune only the parameters a 10 , a 20 and a 01 to fit the hadron observables at vanishing temperature and magnetic field. The fitted results of the hadron properties are listed in Table II. In the fitting process, the f π and f K are obtained via the minimums of Γ k=0 . Then the constituent masses of The masses of the mesons are given by the eigenvalues of the Hessian matrix H i,j The specific form of the Hessian matrix is quite complicated and can be found in Refs. [68,88]. The Polyakov-loop potential V poly also needs to be specified for solving Eq. (16). We use the polynomial form to parametrize V poly [101], which reads where t = (T − T YM )/T YM . T YM is the critical temperature for the deconfinement phase transition in pure Yang-Mills theory. When quarks are added, t YM needs to be rescaled to account for the unquenched effect where α and T QCD are obtained by fitting lattice QCD data. In this work we take the parameters recommended in Ref. [93] which are collected in Table III. Note that we are working at zero chemical potential, the Φ remains thus real. We set thenΦ = Φ in the rest of this article. There is one more thing we should mention. When the external parameters (such as 2πT or eB) are comparable to the cutoff Λ ∼ 1 GeV, the initial effective action Γ k=Λ will grow a dependence on these external parameters. This dependence is omitted in the above discussion and we employ the method proposed in Refs. [92,93] to compensate for this neglected dependence. More detailed discussions can be found in Ref. [102]. The basic idea is to start the flow from a high enough scale in which all the external parameter dependence can be neglected, and then flow down to the original cutoff Λ. As a result, the initial condition Γ Λ will get modified where s strands for all the external parameters (T, eB in this article). ∆Γ Λ [s] represents the external parameter effect and reads The second line follows from that external parameter dependence can be neglected when the ultraviolet cutoff is large enough, i.e., Γ ∞ (s) = Γ ∞ (0). Note that mesons will become too heavy above Λ and decouple from the system, we only substitute the fermion part of Eq. (12) into Eq. (22). For numerical convenience, we set the temperature to zero in the ∂ k Γ k (0) term of Eq. (22) and keep the magnetic field dependence. The specific form of Eq. (22) now reads where E f = k 2 + m 2 f . Generally, quark mass m f depends on ρ 1 and ρ 2 . However it is shown in Ref. [92] that setting m f to its constituent mass value is a good approximation. Thus we set m u,d = 300 MeV and m s = 430 MeV in our calculations.

IV. CALCULATIONS AND RESULTS
Once the effective action Γ is known, many thermodynamic quantities of the system can be obtained. Some quantities are easy to get while some others contain cutoff artifacts and should be handled carefully. We will describe how to obtain these quantities below.
First of all, the minimums of the Γ(ρ 1 , ρ 2 , Φ) are related to the order parameters and we useρ 1 ,ρ 2 ,Φ to denote them. We takeρ 1 as the order parameter for chiral phase transition becauseρ 1 is related to the condensates of light quark and strange quark, which reads And it is usual to useΦ as the order parameter for deconfinement phase transition. We define the pseudo-critical temperature as the inflection point of these order parameters where they change most rapidly.
The pressure P and the entropy density S can be derived directly from Γ as Meanwhile, we can expand the pressure P (T, eB) in terms of eB as whereχ(T ) is the bare magnetic susceptibility which is usually taken to characterize the leading response O((eB) 2 ) of the pressure to the magnetic field. Unfortunately, both P (T, 0) and bare magnetic susceptibilitỹ χ(T ) contain cutoff artifacts and need to be renormalized. We follow the procedures taken in Refs. [27][28][29] to define the renormalized magnetic susceptibility as and χ(T ) is obtained by a polynomial fitting to ∆P [68] in the weak magnetic field region (eB < 0.1 GeV 2 ), which reads Cutoff artifacts are all removed by the subtractions in ∆P .
Magnetic susceptibility χ(T ) has also been studied by combining the FRG approach with the non-interacting gas approximation in Refs. [68,103], in which the whole system is approximated as non-interacting quark and (pseudo-)scalar meson gases [30]. χ(T ) can then be separated into two parts, namely where χ q (T ) and χ m (T ) refer to the contributions from quarks and mesons respectively. After some calculations, we have the two terms explicitly as where the summation for bosons is over charged mesons (π + , π − , K + , K − , a + 0 , a − 0 , κ + , κ − ). Note that the masses of the particles have a dependence on the temperature and neglecting this dependence will lead to a negative definite magnetic susceptibility at any temperature. It is proposed in Ref. [68] to use the temperaturedependent masses and then include additional logarithm terms in Eq. (30), which are vanishing in the original calculation. These temperature-dependent masses m f (b) (T ) and Polyakov-loop Φ(T ) can be obtained from FRG approach, we can thus substitute the results into the above formulas and do the integrals numerically.
To check our calculations further, Adler function [104,105] is also calculated. It's mentioned in Ref. [27] that χ(T ) is closely related to the Adler function which can be calculated in perturbation theory. The relation can be formulated as where D(2πT ) is the Adler function and the momentum Q of virtual photon is replaced by Matsubara frequency 2πT . This relation enables us to compare our results with those given in perturbation theory. In the following subsections, we will represent and discuss our main results.

A. ORDER PARAMETERS AND PHASE DIAGRAM
The calculated results of the constituent quark masses as functions of temperature at several values of magnetic field strength for the light quark and strange quark are illustrated in Fig. 1. We show that, at a fixed temperature of any value, the constituent quark mass increases with the ascending of the magnetic field strength. This phenomenon is just the so called magnetic catalysis and is explained by the dimension reduction of quarks under external magnetic field, see, e.g., Ref. [47]. In order to take the condensate of light quark and strange quark both into consideration, we use the ρ 1 as the order parameter for chiral phase transition and the obtained variation behavior against the temperature is illustrated in Fig. 2. The T χ c extracted from the inflection point of the ρ 1 at eB = 0 is T χ c = 199 MeV. Figs. 1 and 2 both show evidently that the 2+1 flavor PQM model with FRG approach still displays the character of magnetic catalysis and more sophisticated models are needed to realize the inverse magnetic catalysis.
The restoration of chiral symmetry can also be seen from the temperature dependence of the mass spectrum of mesons displayed in Fig. 3. We can see evidently that the masses of π and σ mesons get close as the temperature increases and become finally degenerate, which signals the restoration of chiral symmetry. It is also apparent that the mesons become much heavier than two quarks when chiral symmetry is restored, this means that the mesons become unphysical degrees of freedom as the chiral symmetry is restored since the mesons have in fact melted to quarks and antiquarks.  As usual, we take the Polyakov-loop Φ as the order parameter for deconfinement phase transition and the corresponding pseudo-critical temperature T d c is defined by the inflection point of the Φ. The obtained temperature dependence of the Φ at several values of magnetic field strength is shown in Fig. 4. The obtained T d c in case of vanishing magnetic field is 182 MeV. And the T d c decreases slightly as the magnetic field strength increases, which is consistent with the result given in Ref. [59]. Such a slightly decreasing behavior of the T d c can be explained by the energy change of quarks. When the magnetic field is turned on, the dispersion relation of the quarks get modified where l labels the Landau level and s = ±1 denotes the spin orientation of the particle. Note that quark can stay at the lowest Landau level (l = 0 and s = 1) and its energy is reduced, that's exactly how dimension reduc- tion happens. It will take less energy to create a quark from vacuum, then the Polyakov-loop Φ increases since Polyakov-loop Φ ∼ e −βFq where F q is the free energy of a single-quark system [106]. And in turn, the T d c gets lowered. We will discuss the response of different particles to the magnetic field further in the following subsections. Here, we would like to mention that, since the magnetic field has not yet been included explicitly in the Polyakovloop potential V poly in this work, the external magnetic field affects the Polyakov-loop rather slightly, which is consistent with what has been given in Ref. [59].
In Fig. 5, we display the obtained magnetic field strength dependence of the T d c and T χ c . It indicates clearly that the T χ c and T d c split from each other as the magnetic field increases. Considering the increasing amplitude of the T χ c with respect to the magnetic field strength, our present result is almost exactly the same as the previous result in FRG approach [67] in which T χ increases about 40 MeV at eB = 20 m 2 π . And in meanfield approximation (e.g., in QM model), the increasing amplitude is a little mild, Ref. [87] gives about 30 MeV increasement at eB = 20 m 2 π .

B. PRESSURE AND ENTROPY DENSITY
In Fig. 6, we show the calculated normalized pressure and normalized entropy density at eB = 0 and compare them with the lattice QCD simulation results [19,20]. In order to make a model-independent comparison, the temperature axis is normalized by the T χ c calculated in lattice QCD (T χ c = 154 MeV [19]) and the present PQM + FRG approach (T χ c = 199 MeV), respectively. We can see from the figure that our results are in reasonable agreement with lattice QCD results after the temperature normalization is made. The result shows evidently that the entropy density increases generally with the rising of temperature. Furthermore, the monotonicity of the normalized entropy density function in terms of the temperature changes from up-concave to up-convex as the temperature increases to a specific value. With the chiral pseudo-critical temperature fixed by analyzing the quark condensate (and the constituent quark mass) in mind, one can notice that such a inflection temperature for the normalized entropy density is almost the same as the T χ c for the chiral symmetry to be restored. When the magnetic field is turned on, many thermodynamic quantities (such as pressure) will develop a dependence on the eB and consequently their normalized values will diverge at low temperature [27]. However, entropy density S doesn't contain such a divergence and S need not to be renormalized. We will thus analyze the normalized magnetic field induced entropy density [S(T, eB) − S(T, 0)]/T 3 in the following.
The calculated variation behaviors of the normalized magnetic field induced entropy density in terms of the temperature at several values of the magnetic field  [19,20]. strength are displayed in Fig. 7. One can recognize easily from Fig. 7 that the entropy density decreases with the increasing of the magnetic field strength in the temperature region in which the chiral symmetry is dynamically broken. Nevertheless, it increases as the magnetic field strength gets stronger when the chiral symmetry is restored in high temperature region, which is qualitatively consistent with the results of lattice QCD [27] and HRG model [30]. This behavior can be explained as below. From Figs. 1 and 2, we can see that the magnetic field strengthens the quark condensate drastically at low temperature. It implies that the QCD system at low temperature selects out the microscopic states whose condensate doesn't vanish when the magnetic field gets stronger. This selection makes the system more ordered and, in turn, the entropy density decreases. At high temperature, this selection becomes less important and other effects should be considered to account for the increase of the entropy density. More precisely speaking, if we approximate the system as non-interacting gases, the entropy of the PQM system receives contributions from both quarks and mesons. Combining Eqs. (25) and (30), we can see that the entropy density of the charged (pseudo-)scalar mesons (such as π ± , · · · ) and free quarks response differently to the magnetic field: the energy of each π ± meson is increased by the magnetic field (see Eq. (32)), thus there are less microscopic states available to the π ± system at certain temperature, which leads to a lower entropy density; while the energy of each quark is decreased, hence the entropy of quark system increases. Since the induced entropy of the whole system is mainly determined by the lightest charged mesons π ± when the chiral symmetry is broken, the entropy density decreases then at low temperature as the magnetic field gets stronger. When the temperature is high enough, chiral symmetry is restored and mesons get melted. The entropy density of the system is determined mainly by quarks and thus increases with the increasing of the magnetic field. We will try to discuss these behaviors from a more fundamental point of view in Sec. IV C.
Moreover, we should note that the calculated entropy density at high temperature ( > ∼ 120 MeV) still shows quantitative difference with that given in lattice QCD simulations [27]. Our presently calculated normalized magnetic field induced entropy density at different magnetic field strengths are quite small compared with the results given in Ref. [27]. This may be due to the missing of the inverse magnetic catalysis effect in PQM model. More sophisticated model is required to solve this discrepancy.

C. MAGNETIC SUSCEPTIBILITY AND ADLER FUNCTION
We collect the χs obtained through the present (P)QM+FRG approach, non-interacting gas approxima-tion (see Eq. (30)) and lattice QCD simulations together in Fig. 8. The QM+FRG approach here also takes the parameters in Table I except that the Polyakov-loop and its potential are removed. The non-interacting gas approach approximates the whole system as non-interacting quark and (pseudo-)scalar meson gas, and the particles contribute individually to the magnetic susceptibility [68] as shown in Eq. (30). Moreover, the temperature-dependent particle masses and Polyakov-loop in Eq. (30) are taken from the present PQM+FRG results (as displayed in Figs. 3 and 4). The figure manifests apparently that, at temperature below T χ c , χ is negative and QCD system is diamagnetic. This is because χ receives contributions mainly from the lightest pseudo-scalar mesons π ± below T χ c . It's well known that spin-less charged boson displays diamagnetic character due to its orbital motion perpendicular to the magnetic field. When the temperature exceeds T χ c , mesons become too heavy and decouple from the system, in fact, get melted to quarks. χ will be solely determined by quarks then. External magnetic field will align quark's magnetic moment parallel to it and therefore leads to paramagnetic at high temperature. We can see from the above discussions that the transition of the magnetism of QCD system is closely related to the chiral symmetry, therefore the transition temperature of magnetism is near T χ c . Considering Eq. (28) and the consistency between our present result and those given in lattice QCD simulations [27][28][29], one can recognize that our presently obtained magnetic field induced pressure also agrees with the lattice QCD result [29] very well. Comparing the results via the FRG approach embodied with the PQM and the QM models, we can see that PQM model predicts a higher χ than QM model in the temperature region above the T χ c and a much lower χ below the T χ c . This is because of the effect of Polyakovloop. Quarks are confined at low temperature in PQM model and their contributions to χ are suppressed, which leads to a lower χ. At high temperature, gluon's effects are taken into consideration in PQM model through the Polyakov potential V poly , which can account for the higher χ. The difference between the PQM and the QM models indicates that the effect of Polyakov-loop is important for calculating the magnetic susceptibility quantitatively. We also note that the result from noninteracting gas approximation seems agree with lattice QCD simulations quite well, but the increasing rate is larger than those given in lattice QCD at the temperature higher than the T χ c . To check our calculation results and compare with the lattice QCD simulation result further, we show the calculated Adler function and the comparison with the perturbative QCD result [105] and lattice QCD result [27] in Fig. 9. The figure shows clearly that, besides the reasonable agreement between our presently calculated result and the perturbative QCD result in high temperature region, the drastic variation behavior around the pseudocritical temperature also coincides with the lattice QCD result. We also note that the calculated Adler functions from PQM and QM model behave differently. QM model only consists of free quarks at high temperature, thus its Adler function is close to the perturbative result at the lowest order O(α 0 s ). While quarks in PQM model are still coupled with gluons at high temperature by the effect of Polyakov-loop, thus the corresponding Adler function is more close to the higher order result at high temperature. Before we close this section, we try to explain the re-sponse of QCD system to the external magnetic field in an intuitive picture. At low temperature, chiral symmetry is dynamically broken and quarks prefer to come into pairs, then the whole system is composed of various quark pairs (in fact, quark-antiquark pairs, here). When the magnetic field is turned on, these quark pairs will rotate in the plane perpendicular to the magnetic field in a way that produces anti-parallel magnetic moment to the external field, which leads to a negative magnetic susceptibility. Interestingly, these behaviors of quarks have many similarities with the electrons in BCS theory [107], in which electrons form Cooper pairs and lead to Meissner effect. Moreover, these quark pairs can be equivalently regarded as various mesons and we can see from the previous discussions that the entropy of these quark pairs decreases with the increasing of the magnetic field. At high temperature, the chiral symmetry gets restored. Quarks stop forming pairs and the whole system contains only quarks. Then, external magnetic field will align the magnetic moment of the quarks parallel to it and leads to a positive magnetic susceptibility. The increasing feature of the the positive magnetic susceptibility with respect to the temperature means that there exist more magnetic moments need to be aligned. As a consequence, the entropy density of the quark system increases with the increasing of the magnetic field strength. These characteristics can be regarded as an intuitive picture for the magnetic catalysis in simple phenomenological models.

V. SUMMARY AND REMARKS
In this paper, we investigated the thermodynamical properties of the strong interaction matter in an external magnetic field by applying the FRG approach to (2+1) flavor PQM model. Our calculation results of the temperature dependence of the pressure and the entropy density are qualitatively consistent with the lattice QCD result while quantitative discrepancy still exists when the magnetic field is turned on. The calculated constituent masses of light quark and strange quark indicate that the pseudo-critical temperature for the chiral symmetry to be restored increases with the increasing of the magnetic field strength, which is consistent with the magnetic catalysis effect.
We have also studied the response of the QCD system to the magnetic field by analyzing the magnetic field induced entropy density, the magnetic susceptibility and Adler function. Our calculated results indicate that the property of QCD system undergoes a transition with the increasing of temperature: at temperature below T χ c , the magnetic field reduces the entropy density and the pressure of the system; after chiral symmetry is restored, the magnetic field begins to increase the entropy density and the pressure. This transition is closely related to the chiral symmetry breaking and its restoration. Especially, an intuitive picture for the magnetic catalysis in simple phenomenological models is given.
In addition, because the model we have taken is simple so that the gauge field effects and the feedback effects of the magnetic field on the gauge interaction (see, e.g., Ref. [52]) have not yet been taken into account, our present result is then valid only to the pure and direct effect of the magnetic field on the strong interaction matter, but could not rule out the possibility of the inverse magnetic catalysis effect, at least the retarded magnetic catalysis [62,71,78]. In other words, our investigation indicates that the underlying physics responsible for the inverse magnetic catalysis (or magnetic inhibition) is beyond the physics demonstrated in the PQM model and other similar phenomenological models. To include the gauge field effects and explore the physics for the magnetic inhibition, one may take the nonlocal NJL model [108], or the Dyson-Schwinger equations of QCD [78], or the approaches taking the feedback effect of the magnetic field on the strong interaction to carry out the calculations. The related investigation is under progress.