The Evolution of Meson Masses in a Strong Magnetic Field

Spectra of $q \bar{q}$ hadrons are investigated in the framework of the Hamiltonian obtained from the relativistic path integral in external homogeneous magnetic field. The spectra of all 12 spin-isospin s-wave states, generated by $\pi$ and $\rho$ mesons with different spin projections, are studied both analytically and numerically on the lattice as functions of (magnetic field) $eB$. Results are in agreement and demonstrate three types of behavior, with characteristic splittings predicted by the theory.


Introduction
The influence of magnetic field (MF) on the strong interacting particles is an actively discussed topic, see, e.g. a recent review [1]. When MF is not ultra-intense (eB ≪ σ, where σ = 0.18 GeV 2 is a confinement string tension) 1 , the main characteristics related to the behavior of hadrons in MF are magnetic moments and magnetic susceptibilities, while for the strong MF limit (eB ≥ σ) the hadron energy and width depend on MF directly.
The PIH method has appeared to be well suited to the inclusion of an arbitrary external MF. Here one obtains simple expressions for magnetic moments of hadrons, mesons [52] and baryons [53], which are in a good agreement with available experimental and lattice data, as well as with existing model calculations. We stress at this point, that in all calculations done within the PIH framework, the final results are expressed in terms of basic QCD parameters -string tension σ, α s and current quark masses m q .
A sample of light neutral meson masses in MF (actually, the meson energies for zero longitudial momentum) has been calculated with PIH framework in [40][41][42], and the Nambu-Goldstone (NG) modes in MF have been studied in [47]. In all cases the resulting values of M i (B) are in reasonable agreement with lattice data from [13]- [14]. The three-body neutral systems in strong MF were studied with PIH in [54], but there is no lattice data now to compare with.
In a general case, solving the spectral problem for hadrons in MF is a cumbersome task. To proceed with analytic calculations, one should use some special techniques. One of them is the Pseudomomentum approach. It was introduced in [55] to separate center-of-mass (c.m.) motion from the relative motion in the nonrelativistic Hamiltonian for the neutral system in MF. This approach was extended to the relativistic sector in the PIH framework for two-body systems in [40][41][42] and for three-body systems in [56]. The Pseudomomentum approach is applicable only for electrically neutral systems, and for the charged ones an exact analytical answer was obtained only in an unphysical model of charged meson with equally charged quark constituents [41].
Below we are suggesting a new approximate analytic method of constituent separation (CS) that allows to get a quantitative result for any meson masses with 15% accuracy for the strong MF (eB ≫ σ) and with 20% accuracy for eB < σ. As will be shown, the CS method allows to study charged and neutral systems in the same way. To introduce it, we first write the relativistic Hamiltonian in MF within PIH formalism and exploit the oscillator representation for the confinement interaction used before in [40][41][42] with 5% accuracy. This allows to split the Hamiltonian into transversal and longitudial (with respect to the MF direction) parts analytically. All the rest interaction -one-gluon exchange, spin-dependent and self-energy interactions are studied perturbatively.
Our final results for the neutral mesons in MF are obtained in two independent ways: via Pseudomomentum and the CS methods, which allows to check the accuracy of our results.
The paper is organized as follows: in Section 2 we write relativistic Hamiltonian and discuss the main features of CS method. (Details of this method are discussed in Appendix A). As a result we obtain in Section 2 the hadron mass and the ground state wave function as a function of eB and σ for an arbitrary meson in MF. In Section 3 a classification of meson mass trajectories with different spin and isospin projections is given with the corresponding asymptotics in high MF regime eB → ∞. In Section 4 the perturbative correction due to the one-gluon exchange is calculated and the absence of the color Coulomb collapse is demonstrated. The CS wave function for neutral mesons is discussed in in the Appendix B. In Section 5 the spin-spin interaction in MF and the seemingly possible "hyperfine collapse" is discussed. In Section 6 a general discussion of the we spin-isospin splitting is given. In Section 7 we study the chiral and nonchiral treatment of pion masses in MF. In Section 8 the details of our lattice calculations are given. Results of both analytic and lattice results are discussed in the concluding Section 9.

The relativistic Hamiltonian of quark systems
We start from the relativistic Hamiltonian of the N-quark system in an external homogeneous MF, which according to [40]- [45] is where ω i are virtual quark energies to be intergrated over in the path integral, and m q i are current quark masses. At this step we neglect any internal interactions between quarks, i.e. confinement, gluon-exchange, etc. It is convenient to choose symmetrical gauge for MF A i = 1 2 (B × r i ) which allows to define an angular momentum projection m i for each quark as a quantum number. The spectrum of (2.1) with m i = 0 is According to [40]- [45] the physical spectrum is given by the stationary point value of ε, with respect to ω i It is easy to see that this spectrum coincides with the solution of the Dirac equation for N non-interacting relativistic particles in MF. As in [41] we now introduce the confining interaction V conf , which is treated nonperturbatively, while the other interactions like one-gluon exchange V OGE , spin-dependent interaction a SS and self-energy corrections ∆M SE are treated perturbatively in the next sections. The Hamiltonian becomes with the ground state eigenvalue M d (nonperturbative, or dynamical mass) and the ground state wave function |Ψ 0 . The total meson mass is a sum of M d and the perturbative corrections One can note, that the contribution of the V conf in strong MF (eB ≫ σ) is negligible in the plane transverse to the MF direction and should be retained only for lowest levels, which we call "zero hadron states" (ZHS) (see below). Another feature is that in strong MF regime the translational invariance of the center-of-mass (c.m.) is broken due to magnetic forces (each quarks is placed on its own Landau level), but the confinement still defines the motion of quarks in the direction along the MF.
To simplify calculations we chose the confining term in the variable quadratic form [41,42], restoring its original linear form at the stationary point (it was checked to be accurate within about 5%), namely where γ is variational parameter and σ = 0.18 GeV 2 is a confinement string tension. The dependence of the string tension σ on the MF is caused by the fluctuating qq pairs embedded to the string and provides a correction about ∆σ σ ∼ 15% at eB ∼ 1 GeV 2 . This phenomenon was studied on the lattice in [43] and was confirmed within PIH formalism in [44]. The correction to the ground state caused by this effect is beyond the declared accuracy and is neglected in what follows. To produce an approximation for the energy, one should minimize the resulting state energy obtained from the Hamiltonian (2.4) with respect to ω i and γ simultaneously.
The oscillator approximation (2.6) gives an advantage to separate motion along the z axis(parallel to the MF) and in x − y plane where the motion along the z-axis is defined by the Hamiltonian where we use c.m. reference frame with P 3 = p The longitudial part of the ground state energy is For the motion in the transversal plane one can use an approximation of decoupled quarks at large MF, making the following substitution (r where c.m. position r 0 ⊥ is fixed at the origin in x − y plane. This approximation corresponds to the configuration where the confinig string connects each quark to the c.m., i.e. the string is effectively elongated. The magnetic energy of each quark in strong MF (Landau level) is larger than the confinig interaction with the factor eB σ , which make this approximation legitimate at eB > σ regime. To extend our method to the eB < σ region, where the behaviour is mostly defined by confinement, one should introduce an effective sting tension σ 1 and σ 2 for each part of the string, connecting quarks to the c.m.
to compensate an effective string elongation. As shown in Appendix A, the appropriate values of σ 1 , σ 2 are (2.12) Using this "σ-renormalization" procedure, one can show that the dynamical mass of the ground state M d , calculated in [41] with the Pseudomomentum technique for neutral mesons, exactly coincides with the dynamical mass obtained in the above CS formalism for the arbitrary value of MF.As a result this approximation make quarks effectively decoupled in x − y plane and one can write The transversal part of the hamiltonian H ⊥ has the ground state energy where σ i are given by (2.12). The total dynamical mass is given by the sum The actual trajectories for the dynamical mass in MF, M d (eB) are obtained using the stationary point conditions in a similar way as (2.3) The corresponding wave function for the ground stateM d is where Ω 1 , Ω 2 and Ω z are harmonic oscillator frequences Comparing (2.17) with the same wave function obtained in [41] for neutral mesons one can see that now we have two elongated ellipsoids for each quark instead of one ellipsoid in η = r 1 − r 2 , but the resulting spectra coincide.

Meson trajectories in MF
We turn now to the general structure of the meson spectrum and the limits of weak (eB < σ) and strong (eB ≫ σ) MF. For small MF both γ and ω are independent of MF in the leading order and the lowest order the correction to the dynamical mass is where µ is the magnetic moment of the hadron, and the c|eB| term is c.m. energy contribution (the lowest Landau level) in MF for the charged mesons (note, that in this paper we discuss only s-wave hadrons and all orbital momenta are zero). Magnetic moments in PIH formalism have been calculated in [52] for mesons and are in good agreement with experiment and available lattice data. It is easy to see there that for massless quarks the expansion in (3.1) is actually done in powers of eB σ . For the strong MF limit the situation is more complicated. Confining ourselves to the lowest Landau levels (LLL) for all quarks and antiquarks, i.e. n (i) ⊥ = 0 in (2.9) and (2.14), we can separate out the hadrons, which consist of only LLL states of both quarks with e i σ i z = |e i |, i = 1, 2. These states are MF-independent at eB → ∞ and we shall call them "zero hadron states" (ZHS). Note, that ZHS do not possess definite total spin and isospin quantum numbers.
All other hadron states, except for ZHS, will have energies growing with MF as ∼ |eB| and therefore thermodynamically suppressed at large MF. In the limit of strong MF the dynamical masses for ZHS can be written as The stationary point analysis according to (2.16) for m 1 = m 2 = 0, e i σ i z = 1, i = 1, 2, yields The same result was obtained in [40][41][42]54] with the Pseudomomentum technique.
We turn now to the meson states still with zero orbital momentum and not belonging to the ZHS states, i.e. violating the equality e i |e i | σ i z = 1. The resulting meson energy according to (2.15 The same result occurs when e 1 |e 1 | σ 1 z = −1, e 2 |e 2 | σ 2 z = 1 with replacement e 2 → e 1 Now we turn to the case when both products e 1 |e 1 | σ 1 z and e 2 |e 2 | σ 2 z are negative. In this case one obtains Thus we have three different asymptotic modes for s-wave meson dynamical masses M d in MF, classified with respect to spin projections We shall return to this classification later in Section 6 in our study of spin-isospin splittings in weak MF regime.

One-gluon exchange in MF
The first order perturbation correction for one-gluon exchange potential(OGE, or color Coulomb interaction) in MF entering in (2.5) according to [42] is 3 GeV , and the M 2 B = 2πσ = 1.1 GeV 2 preventing Landau singularity was calculated in [51]. Form (4.1) includes screening of the OGE potential by the quark-antiquark pairs created in MF. This effect prevents the "fall-to-the-center" phenomenon for ZHS hadrons in MF, as shown in the Fig.1 and 7. One can see that the matrix element Ψ|V OGE |Ψ for meson saturates at eB ∼ 10 GeV 2 and the system becomes "asymptotically free" in eB → ∞ limit when Ψ|V OGE |Ψ → 0. The driving force of the Coulomb collapse is an uncontrollable growth of the Coulomb interaction when the system is squeezed by MF forces. The role of screening of the Coulomb interaction in MF has a long story, see e.g. [4][5][6] for atomic systems. The next step is to average the potential (4.1) over the wave function (2.17) obtained by the CS method.
Separating the integration in x − y plane and in z-direction, one has In the momentum space one obtains where Ω i and Ω z are given by (2.18). Comparing this result in case of the neutral meson with the exact one, obtained with Pseudomomentum procedure, one has to make a correction for the wave function, see Appendix B for details.

Spin-dependent corrections
A detailed review of the spin-dependent forces in PIH framework is given in [57]. Here we only emphasize that the spin-dependent perturbative corrections arise from the σ i F · σ j F correlators, where σ i are Clifford 4 × 4 σ µν for i-th quark constituent and F are non-abelian field strength tensors. Averaging over the stochastic gluonic background field, one has two types of corrections -the self-energy term for i = j and color-magnetic spin-spin interaction terms for i = j, where i, j are quark numbers The self-energy correction ∆M SE in (5.1) was used in a large number of calculations [59], confirmed by the experimental data and lattice simulations. In case of an external MF we retain in ∆M SE the value ω A different story is for the spin-spin interaction in (5.1). As it was shown in [7,41,57], the wave function of hadronic and atomic systems becomes "focused" at the origin by MF, i.e. |Ψ 0 (0)| 2 ∼ eB for large MF value. This "magnetic focusing" phenomenon could induce the fall-to-the-center phenomenon for the lowest lying ZHS states. However, as shown in [57], the colormagnetic fields cannot violate the positivity of the qq spectra, implying that some sort of the cut-off parameter must occur in the whole perturbative series with nonperturbative background. Moreover, PIH method has a natural dimensional cutoff parameter for color field λ ≃ 1 GeV −1 -correlation length of the vacuum gluonic background, which should be used to smear δ-function in (5.1) and after the averaging with the CS meson wave function (2.17) one obtains the spin-spin matrix element Smearing procedure prevents the collapse of the meson in strong MF and it stops the unbounded fall of the total mass value in increasing MF. It is important to notice here that the approximation of the confinement potential by the harmonic oscillator potential (2.6) gives too small value for |Ψ 0 (0)| 2 and the hyperfine splitting ∆E = 4 a SS between the non-chiral π − and ρ − mesons at eB = 0 is too small (see Fig. 3 and 7) as compared with realisitc case of linear interaction. Moreover, the pion mass at eB = 0 is additionally shifted down by chiral dynamics, which we shall take into account in Section 7.
6 Spin-isospin splittings in MF As pointed out in Section 1, MF violates spin and isospin symmetries, therefore π 0 , ρ 0 split into 8 states and each π + , ρ + and π − , ρ − states split into 4 states in MF correspondingly. Using the asymptotics (3.7), obtained in Section 3 for strong MF regime, one has 1) ρ + (s z = 1) = |u ↑d ↑ ZHS 2) ρ + (s z = −1) = |u ↓d ↓ II) Here on the l.h.s we have the standard spin-isospin configurations for mesons at zero MF, and on the r.h.s we have asymptotic classification according to (3.7) in strong MF for the corresponding states. The states 1)-4), 5)-8) and 9)-12) are composed of quarks and antiquarks in the combinations which yield the required spin and isospin values of π and ρ mesons at eB = 0.
With increasing MF the eigenvalues of the total Hamiltonian (2.5) at nonzero MF demonstrate two types of phenomena: a) the mixing effect, due to spin-spin forces, equivalent to the Stern-Gerlach phenomenon, when the MF eigenstate can be expanded in two eB = 0 eigenstates; b) the splitting effect, when the zero MF state composed of uū and dd components, splits into two trajectories due to isospin flavor. Finally, the trajectories for charged mesons like ρ + (s z = 1) and ρ + (s z = −1) starting at the same mass at eB = 0, split into two for eB > 0.
To take into account the spin-spin interaction, we choose the basis states | + + , | + − , | − + , | − − in spin space. The states 1) and 2), that corresponds to ρ + (s z = 1) and ρ + (s z = −1) mesons at eB = 0 correspondingly, are diagonal and their dynamical masses are After the stationary point analysis (2.16) one has two sets of parameters (ω , ω ). The total mass of these states according to PIH formalism are given by M total (ρ + (s z ± 1)) gives rise to two trajectories in MF starting at ρ + meson mass at eB = 0. The behavior of states 3) and 4) corresponding to ρ + (s z = 0) and π + at zero MF is more complicated. These states are composed of |u ↓d ↑ = | − + and |u ↑d ↓ = | + − combinations at eB = 0. When the MF increases, the states start to mix in the mutually orthogonal combinations The basis vectors are equal to the π + and ρ + (s z = 0) states at eB = 0. The mixing phenomenon is defined by the non-diagonal spin-spin matrix elements The dynamical masses and the parameters (ω  and to calculate mixing coefficients α, β(ρ + (s z = 0)) = The trajectory E 1 with "+" sign in (6.9) starts from the ρ + (s z = 0) mass at eB = 0 and grows with MF, and the trajectory E 2 with with "-" sign corresponds to the π + at zero MF. The states E 1 and E 2 are mixtures of π + and ρ + (s z = 0) at eB = 0 with mixing coefficients defined by (6.8).

Pion chiral degrees of freedom in MF
Unlike ρ mesons, the pions obey the chiral dynamics and therefore one should take into account how it changes under the influence of MF. This topic was studied in [47] and here we exploit the results of [47] for neutral and charged pions. The most important feature of these results is that the GMOR relations [61] are kept valid for neutral pions in arbitrary strong MF, while they are violated for charged pions. At the same time at zero and small MF (eB < f 2 π ) the pion mass is defined by GMOR relations both in the neutral and charged case, m 2 π = mq| qq | f 2 π . This last dependence m 2 π ∼ O(m q ) defines the main difference between chiral and nonchiral pion trajectories.
We start with the neutral pion case, with the standard GMOR relations.
Since m −+ n,i is fast growing with eB, one can retain in the sums (7.2),(7.3) only the (+−) terms, and obtain as in [47] the asymptotic behavior of the π 0 mass as m 2 π 0 = m M (0) (m (+−) ) 2 , wherem (+−) is close to the lowest mass m n,i . It can be seen that all formulae in derivation of the ECL in [47] are diagonal in isospin flavor, so one can write an independent GMOR relation for each π 0 (uū) and π 0 (dd) mesons, that should split in MF according to (6.1) The result of the calculation for π 0 (uū) (solid line) according to (7.4) is shown in the Fig.2.
In the case of the charged pions, π + and π − , the situation is drastically different, since they loose their chiral properties at large eB > σ, and their asymptotics is defined by the independent u andd quarks (for π + ), the mode I in (3.7), splitted in two trajectories, M +− (eB ≫ σ) ≈ 2 3 eB for π + and M −+ (eB ≫ σ) ≈ 4 3 eB for ρ + due to the π + − ρ + mixing effect mentioned in Section 6. Taking into account the GMOR relations for eB = 0(eB < f 2 π ), the asymptotics for the charged pions can be written as where m 2 π (0) is the pion mass at eB = 0. In Fig.3 we plot the trajectory of M +− (B) and our lattice data together with the lattice data from [28].
One can summarize our method of chiral meson mass calculation as follows. First, one calculates the spectrum of non-chiral ("spectator") meson masses m (+−) n,i in the MF. Second, one uses (7.2), (7.3) to obtain f 2 π and qq as functions of eB. Finally, one exploits GMOR relations to extract the resulting chiral mass dependence on eB. The formalism of (7.2), (7.3) was checked without MF in [63,64], the resulting MF dependence of f 2 π and qq was checked vs. lattice data in [47,65]. Summarizing the results for neutral and charged pions one can say that our theoretical predictions are supported by lattice data, and indeed charged and neutral pions behave quite differently at large MF, violating and not violating respectively the GMOR relations.

Lattice calculations
The ground state energies of π and ρ mesons are calculated in SU(3) lattice gauge theory without dynamical quarks. Technical details were presented in our previous work [19,20]. πlattice a=0.115 πlattice Bali πchiral πnon-chiral Figure 3.
Mass evolution of the chiral (solid line) and nonchiral (dashed line) π − meson in MF in comparison with the lattice data: triangles [28] and quadrangles with lattice spacing a = 0.115 f m (present paper). numerically we found eigenfunctions ψ k and eigenvalues λ k for a quark in the background gauge field A µ . Two types of quarks u and d are considered, which are degenerate in mass. An abelian MF interacts with quarks, so U (1) gauge field is introduced into the lattice version of the Dirac operator D [26] Quark fields obey periodic boundary conditions in space and antiperiodic boundary conditions in time. The MF is quantized in a finite lattice volume. Its value is determined by the following formula where e is the elementary charge. Taking the average over the background field A we introduce the correlators in coordinate space where O 1 , O 2 = γ 5 , γ µ are Dirac matrices, µ, ν = 1, .., 4 are Lorentz indices, x and y are lattice coordinates. We performed the numerical Fourier transform of (8.4) in spatial discrete coordinates and set p = 0 since we are interested in the meson ground state energy. To obtain the masses we expand the correlation functionC(n t ) into the exponential series When the lattice time n t is large, the main contribution to the correlator (8.5) comes from the ground state. Due to the periodic boundary conditions the correlator has the following form where A 0 is a constant, E 0 is the ground state energy, a is the lattice spacing. We find the energy E 0 , as a fit parameter, fitting the lattice correlators by formula (8.6). In order to minimize the errors and to exclude the contribution of excited states we take various values of n t from the interval 5 ≤ n t ≤ N T − 5. The energy of the charged pion is calculated from the correlation function Fig.3 shows the energy of a charged pion for the lattice volume 18 4 , the lattice spacing 0.115 fm and the quark masses 17.13 MeV, which corresponds to the pion mass m π = 395 MeV at zero MF. The resultant energy increases with the MF value. Errors were obtained through the χ 2 method. According to the exponential fall of the correlator (8.6) at large E, the absolute error of E should grow with energy. However, we do not see this tendency clearly since different number of gauge configurations were utilized for different values of the MF. These numbers are shown in Table 1 in case of the charged meson. We observed an increase of the error for the MF values eB > 1.2 GeV 2 , presumably due to the worst convergence of our numerical procedure at high MF values and small quark masses. In Fig.4 we show the correlation functions for different MF values for comparison.
The energy of neutral pion was calculated using the correlation function The resultant energies for the π 0 (uū) and π 0 (dd) configurations diminish with the increase of the MF as shown in Fig.5-6. The correlation functions of charged ρ mesons for three spatial directions are given by the following relations If an abelian MF is directed along the 'z' axis, the ρ meson energy with s z = 0 spin projection to the MF direction is determined by the C V V zz correlator. The ground state energies of the ρ meson with spin projections s z = +1 and s z = −1 are determined by the following combinations of correlators We have obtained that the energy of the ρ − meson with the spin projection s z = −1 diminishes as a function of eB. Fig.7 shows that the energies of ρ − meson with spin projections s z = 0 and s z = +1 increase with the magnetic field value. The energy of the neutral ρ meson was calculated similarly to the charged ρ meson, but in formulae (8.9), (8.10) and (8.11) one has to consider the sum of the correlators for u and d quarks. In

Results and discussion
Our paper contain analytic and numerical lattice results for all 12 π, ρ mass trajectories as functions of eB. The main difficulty with analytic results was threefold. First of all, for charged states there is no universal method of the separation of c.m. and relative coordinates (unlike the case of neutral mesons), and therefore we have used a new special (however approximate, O(15%)) approach, called the CS formalism. Secondly, strong MF in lowest approximation brings in vacuum instability due to OGE forces and due to the hyperfine interaction, both growing fast with eB. We have eliminated the OGE instability taking into account the screening effect [42], as was shown in Section 4, see also Fig.1. For the hyperfine problem we have used the stability theorem of [47], ensuring the nonnegativity of hadron masses in the magnetic field, which implies that in higher order the combined nonperturbative and perturbative effects must stabilize hyperfine interaction. To this end we have employed the smearing radius of the hyperfine term of the order of vacuum correlation length λ ≃ 0.2 f m. Note, that this problem exists also without MF and is usually solved in the same way. Thirdly, the chiral dynamics, which governs pions at zero MF, may be violated by MF, and this was explicitly demonstrated in [52]. Accordingly we had to consider separately charged and neutral pions, where only the latter keep the the chiral properties, see Fig.2 and 3. . Mass evolution of (π 0 , ρ 0 )(uū) quartet in MF from analytic (PIH) and lattice data (black circles are from [27]).
Indeed we show in Fig.5 and 6 the behaviour of the π 0 (uū) and π 0 (dd) masses in MF, which follows from the GMOR relation, where both f π (eB), and qq (eB) are calculated via non-chiral qq eigenvalues in MF. Note, that the chiral π 0 mass in Fig.2 and the nonchiral π 0 mass are similar in behavior but differ in scale. The latter is due to fact, that in chiral dynamics m 2 π is proportional to the quark mass m q , see Eq.(7.4). Moreover, the nonchiral neutral pion mass becomes negative for eB > 0.6 GeV, when the standard hyperfine cut-off of λ ≈ 1 GeV −1 is used, which might require a smaller λ. This fact calls for an additional investigation.
As for charged pions π, one can see in Fig.3 a drastically different behavior which has growing asymptotics of the type I according to (3.7) for chiral (m 2 π = O(m 2 q ), lower curve) and nonchiral (m 2 π = O(σ) higher curve), cases. One can see in Fig.3 a reasonable agreement of lower curve with the lattice data of [14] while our present lattice data in Fig.3 correspond to much larger m q and therefore are shifted upwards.
Turning to the ρ mass trajectories, one must remember our classification in Section 6, which implies, that both π 0 , ρ 0 lines split into (uū) and (dd) species and the growing trajectories are proportional to |e q |B, yielding for those a ratio equal to √ 2. The ρ 0 (s z = 1, −1) PIH trajectories in Fig.5 and 6 agree well with our lattice data and with lattice data from [27], as well as π 0 (uū) and π 0 (dd) trajectories.
Note the difficulty in lattice evaluation of the ρ 0 (s z = 0) lines which mix with the much lower π 0 trajectories.
A very interesting situation occurs for ρ − (ρ + ) mass trajectories, presented in Fig.7. Only one of those belong to the ZHS type and tends to a constant at large MF, and both lattice and analytic curve agree within our accuracy O(15%), approximately the same kind of agreement is seen in Fig.7 for the trajectory of the type I, ρ − (s z = 0) and that of the type II, ρ − (s z = 1). Summarizing, one can conclude, that our lattice data agree with analytic predictions within our accuracy limits and our classification and theory based on the PIH formalism for all s-wave π, ρ mesons give a realistic physical picture in this section.

A String tension renormalization in CS method
The string tension renormalization procedure could be illustrated by an analogy from classical mechanics, where two point masses are connected by the spring with the following classical Lagrangian In what follows we take substitution m 1 = ω 1 , m 2 = ω 2 and k = σ 2γ , which makes the Lagrangian (A.1) canonically conjugated to the relativistic Hamiltonian (2.1) with the confinement potential was taken in oscillator form (2.6) at B = 0 up to momenta-and coordinate-independent terms. The Lagrangian could be canonically quantized in c.m. reference frame where E 0 is ground state. On the other hand, one can describe the same system as two independent oscillators with opposite phases (for P = 0), each of them is connected to the c.m. with its own spring, with Lagrangian where the stiffnesses are k 1 k 2 = m 1 m 2 . To proceed further with canonical quantization procedure, one has for ground the state energy for (A.3) The expression (A.4) should take into account that in the quantum case the phases of two harmonic oscillators should be entangled due to constraint m 1 x 1 + m 2 x 2 = 0. An explicit derivation of the "σ-renormalization" requires Dirac quantization formalism for constrained systems. Here we use an heuristic way, based on the correspondence principle. One can substitute the constraint m 1 x 1 + m 2 x 2 = 0 and The ground state energy (A.7) equals to (A.2) if one redefines k 1 and k 2 as (A.8) The resulting ground state energies, obtained with CS method and "σ-renormalization" procedure exactly coincide with the corresponding energies for the neutral mesons in the Pseudomomentum technique [41] in a whole range of MF, see e.g. (2.12-2.14) in the main text of the paper.

B The wave function of the neutral meson in CS method
One can suppose that the results for the averaged operator (4.3) for neutral mesons coincides with the exact one from [41] with Pseudomomentum technique in strong MF regime (eB ≫ σ) because of one-to-one correspondence for dynamical masses in Section 2. However, there is a difference about 30% (especially for ZHS states) because of the c.m. fixing procedure, i.e. the translational invariance breaking in CS method (see Section 2). The nature of this discrepancy is the lowest Landau level (LLL) degeneracy in angular momentum projection m in symmetric gauge, when MF is strong enough to make the confinig force negligible. Let's consider the Hamiltonian for a single particle in MF to illustrate this statement (B.1) The corresponding spectrum for this Hamiltonian is E = Ω(2n + |m| − m + 1), (B.2) Comparing transversal parts of the wave functions for CS method (2.18) as (f) and Pseudomomentum method as (p) in strong MF for ZHS meson

(B.7)
It's evident that the term exp − eB 4 (x 1 · x 2 ) in (B.7) is formed by the power series in m 1 and m 2 entering before the exponent in (B.4). So, the difference between the CS wave function ψ f and the exact wave function ψ p is given by the superposition of the degenerate LLL basis wave functions. This additional term in (B.7) gives about 30% of the total value of the CS Coulomb integral (4.3) for neutral mesons. The nature of this underestimation is clear -an additional term in (B.4) recovers translational invariance of the c.m. for the ψ f wave function. Also one should note that this correction doesn't exist for the charged meson case because of lack of the c.m. tranlational invariance due to c.m. precession in MF. The final step is to add an additional multiplier exp − eB 4 (x 1 · x 2 ) to our CS wave function (2.18) by hand according to previous speculations. This modification gives us 10% accuracy for the Coulomb correction integral (4.3) in comparison with one was obtained in [41].