Real-time dynamics of Chern-Simons fluctuations near a critical point

The real-time Chern-Simons diffusion rate is studied in $(1+1)$-dimensional massive Schwinger model with a $\theta$-term. We evaluate the real-time correlation function of electric field that represents the topological Chern-Pontryagin number density in $(1+1)$ dimensions. Near the parity-breaking critical point located at $\theta=\pi$ and fermion mass $m$ to coupling $g$ ratio of $m/g \approx 0.33$, we observe a sharp maximum in the Chern-Simons diffusion rate. We interpret this maximum in terms of the growth of critical fluctuations near the critical point, and draw analogies between the massive Schwinger model, QCD near the critical point, and ferroelectrics near the Curie point.

Recently, quantum algorithms have emerged as an efficient (and potentially superior) way to explore the dynamics of quantum field theories, including the Schwinger model . Previously, we have addressed the realtime dynamics of vector current [47] (a (1 + 1) analog of the chiral magnetic current [56,57]) induced by the "chiral quench" -an abrupt change of the θ-angle, or the chiral chemical potential. In this paper, we will explore the connection between the real-time topological fluctuations and criticality using Schwinger model as a testing ground. The massive Schwinger model possesses a quantum phase transition at θ = π between the phases with opposite orientations of the electric field, see Fig.1. For m g, this phase transition is first order. However, it was shown by Coleman [58] that the line of the first order phase transition terminates at some critical value m * , where the phase transition is second order. The position of this critical point was established at m * ≈ 0.33g [59][60][61]; the resulting phase diagram is shown in Fig. 1. The phase diagram of the theory in the (m/g, θ) plane is thus reminiscent of the phase diagram of QCD in the (T, µ) plane of temperature T and baryon chemical potential µ [62][63][64].
To understand the physics behind this phase diagram, we need to recall the role of θ-angle in the model. The action of the massive Schwinger model with θ term in * kazuki7131ATgmail.com † dmitri.kharzeevATstonybrook.edu ‡ ykikuchiATbnl.gov with / D = γ µ (∂ µ − igA µ ). Note that the gauge field A µ and the coupling constant g have mass dimensions 0 and 1, respectively. Upon a chiral transformation, ψ → e iγ5θ ψ andψ →ψe iγ5θ , the action is transformed to, 1 It is clear from (2) that the massive theory with a positive mass m > 0 at θ = π is equivalent to the theory at θ = 0 but with a negative mass −m. This line terminates at m * ≈ 0.33g at the critical point, where the phase transition is second order. For small masses m m * , the electric field is screened by the production of light fermion-antifermion pairs.

II. TOPOLOGICAL FLUCTUATIONS NEAR THE CRITICAL POINT
Let us now discuss topological fluctuations in massive Schwinger model. For that purpose, it will be convenient to use the Hamiltonian formalism with the temporal gauge, A 0 = 0. From the action (2) the canonical momentum conjugate to A 1 can be read off as Π =Ȧ 1 . The corresponding Hamiltonian is then given by In the bosonized description of the theory, the Hamiltonian is given by (see Appendix A for details) The potential of the model is thus given by where µ = g/ √ π is the mass of the scalar boson, and the dimensionless coefficient c is given by c = e γ /2π with the Euler constant γ = 0.5774. The electric field is related to the boson field as where the first term is the quantum, dynamical contribution and the second one is the classical background induced by the θ angle, E cl = −gθ/2π. The physics is periodic as a function of θ in the absence of boundaries -as θ increases, the electric field becomes capable of producing a fermion-antifermion pair (or kinkantikink pair, in bosonic description) that screens it. At θ = π the potential takes the form When m µ, the potential has two well separated minima at ϕ ≈ 0, − √ π, associated with spontaneous symmetry breaking. As follows from (6), the minima ϕ = 0, − √ π correspond to the electric fields E = −g/2 and E = g/2, respectively.
When m µ, there is a single minimum at ϕ = − √ π/2, which according to (6) corresponds to the phase with no electric field. This is because at small m, the electric field is easily screened by the production of fermionantifermion pairs -so we are dealing with the screened phase. At some critical value m ≈ m * , the effective potential becomes flat -this corresponds to the critical point with a second order phase transition.
For m > m * , the minima are separated by a potential barrier, and the transition between them (corresponding to the change in the direction of the electric dipole moment of the system) is first order -for example, having a domain with an opposite orientation of the electric dipole moment would cost an additional energy due to the "surface tension". Due to this extra energy, the fluctuations of the electric dipole moment at m m * are suppressed. At the critical point, where m ≈ m * , the potential barrier between the two minima disappears, and the "surface tension" of the domains with opposite orientations of the electric dipole moments vanishes. Because of this, the fluctuations of the electric dipole moment near the critical point are strongly enhanced. Here one can draw a useful analogy to the physics of ferroelectrics, where the electric susceptibility exhibits critical behavior near the Curie point [65].
When the ratio m/g becomes very small, the electric field is easily screened by the production of fermionantifermion pairs (or kink-antikink pairs in the bosonized description), and the fluctuations of electric dipole moment again become small. Basing on this qualitative picture (that we will confirm below with a more formal treatment), we expect to see a maximum in the electric susceptibility near the critical point. We will show that this is indeed the case.

III. CHERN-SIMONS DIFFUSION AT ZERO TEMPERATURE
The Chern-Simons diffusion rate is defined as the probability of a process changing the Chern-Simons number to occur per unit volume (length) and per unit time: where Q(x) = ∂ µ K µ is the density of the Chern-Pontryagin number given by the divergence of the Chern-Simons current K µ . In (1 + 1) dimensions, K µ = (g/2π) µν A ν , and Q = (g/4π) µν F µν , where F µν is the electromagnetic field strength tensor. The corresponding change of Chern-Simons number is ∆N CS = d 2 xK 0 = (g/2π) d 2 xA 1 ; the Chern-Pontryagin number is given by the electric field: The diffusion rate (8) is the zero frequency and wavelength limit of the symmetrized Wightman correlator of the topological charge: The fluctuation-dissipation theorem relates the symmetrized Wightman correlator to the imaginary part of the retarded Green's function: where the second equality applies at finite temperature that is usually the origin of diffusion. The retarded Green's function is: where Θ(t) is the Heviside's step function. Hence, to extract the diffusion rate all we need to know is the zero momentum and small frequency limit of the retarded Green's function: The imaginary part of the retarded Green's function is associated with the thermal dissipation and therefore is a T-odd function. As a result the leading term in the small frequency expansion has to be linear in ω, with the coefficient that controls the diffusion rate (12) .
We are interested in the correlation function of the topological charge at zero temperature T = 0. In this case the second equality in (10) cannot be used and we have to use instead the zero temperature limit of (10): and IV. LATTICE SCHWINGER MODEL

A. Lattice Hamiltonian
For the purpose of numerical simulation, we place the theory (3) on a spatial lattice. We introduce the staggered fermion χ n and χ † n [66,67], and lattice gauge field operators L n and U n with an integer n labeling a lattice site; the lattice spacing is a. A two-component Dirac fermion ψ = (ψ 1 , ψ 2 ) T is converted to a staggered fermion ψ 1 (ψ 2 ) → χ n / √ a for odd (even) n. The gauge fields are replaced by e −iagA1 → U n and Π → −gL n , that are placed on a link between nth and (n + 1)st sites. The resulting lattice Hamiltonian is, with the Gauss law constraint: We solve the constraint to eliminate the electric field operators, where we have fixed the boundary electric field, L 0 = 0. By enforcing the relation (17), the states are automatically restricted to the physical ones. We furthermore eliminate the link fields U n by the gauge transformation, with In the present work, we limit θ to be 0 or π to study the critical behavior at θ = π. Upon absorbing cos θ = ±1 in the second line of (15) to the fermion mass m, we arrive at the Hamiltonian, It is noted again that the massive theory with a positive mass m > 0 at θ = π is equivalent to the theory at θ = 0 but with a negative mass −m. The Hamiltonian (20) accords with the latter viewpoint and will be used in what follows.

B. Chern-Simons diffusion rate
We compute the retarded Green's function of electric field operators, Recall that the Chern-Pontyagin number density is proportional to the electric field in 1+1 dimensions (8).
For the purpose of calculating the diffusion rate we take the zero-wavelength limit followed by the zero-frequency limit, where Θ(t) is the Heaviside's step function, V is spatial volume (length) andĒ(t) := dxE(t, x)/V is spatial average of electric field operator. Hence, the diffusion rate takes the form, On a lattice, it is given by withL := 1 Note that the temporal integral is also truncated by T (not to be confused with temperature). We numerically compute the dimensionless diffusion rate, with the dimensionless variablest := (ag 2 /2)t andT := (ag 2 /2)T .

V. RESULTS AND DISCUSSION
We have numerically computed the real-time diffusion rate (25) using a Python package QuSpin [68,69]. Our results are shown in Fig. 2 for g = 1 as a function of the fermion mass m for different values of time T and lattice sizes N . The critical masses m * are extracted by fitting the obtained data points by the following function [60,70]: where a, b, c, and m * are the fit parameters; see Appendix B for details. The extracted values of the critical mass m * (in units of the coupling g) are given in Table I.
The most prominent feature of our results is the existence of a sharp peak in the Chern-Simons diffusion rate near the critical point. The position of the critical point extracted from the Chern-Simons diffusion rate m * /g ≈ −0.33 is in accord with previous studies [59][60][61].T The dependence of the diffusion rate on the lattice size can be understood using the finite-size scaling, with the critical exponents of the transverse Ising model. The Z 2 symmetry of lattice Schwinger model at θ = π and broken parity put it in the same universality class as the 1 + 1 dimensional transverse-field Ising model, which in turn is equivalent to 2D Euclidean classical Ising model. The corresponding finite-size scaling analysis is explained in the Appendix B.
In Fig. 3, we show the ratio of the diffusion rates for different values of timeT . We observe that the characteristic time scale of the fluctuation dramatically increases near the critical point -this is the consequence of the "critical slowing down". At the same time, the magnitude of the fluctuations increases, reflecting the growth of the correlation length, as seen from the peak in Chern-Simons diffusion rate -hence, the fluctuations become slower but larger in magnitude.
The sharp peak in the Chern-Simons diffusion rate near the critical point of the phase diagram may have important implications for the search for the critical point of the QCD phase diagram. In particular, we expect that the CME-induced fluctuations in hadron charge asymmetry [56,71,72] can be significantly enhanced near the critical point.
The similarity to ferroelectrics may allow one to simulate the real-time ferroelectric response near the critical point [65]. Given the importance of ferroelectrics for information storage and processing [73], this may have interesting practical applications.
The Hamiltonian formalism that we use is well suited for a quantum simulation. It would be particularly beneficial for the future studies of higher dimensional or more complicated theories such as QCD, where real-time simulations on classical computers are very challenging.  We review the massive Schwinger model at finite θ (see e.g. [74] for a comprehensive review). Let us start with 1+1D Maxwell theory at finite θ: where E = F 01 = ∂ 0 A 1 =Ȧ 1 is an electric field and the temporal gauge A 0 = 0 is taken. The momentum Π conjugate to A 1 is, obeying [Π(x), A 1 (y)] = iδ(x − y). Hence, the Hamiltonian density takes the form, The Gauss operator must vanish, on physical states. Therefore, Π is fixed once its value is specified, Π| boundary = 0 , at one of the two boundaries.
The boundary value 0 is interpreted as a background electric field. Hereafter, we take 0 = 0 with the background electric field being absorbed into θ.
We now add a dynamical Dirac fermion ψ to have, Then, the Gauss law fixes the electric field in the bulk as where 0 = 0 is taken. We call Π the dynamical electric field because it is induced by fermionic excitations as opposed to E =Ȧ 1 = Π − gθ/(2π), where the second term represents the background electric field.
In line with the main text, we will briefly describe the bosonized form of the theory. The bosonization dictio-nary is as follows: with a dimensionless constant c = e γ /(2π) and the Euler constant γ. Therefore, the action (A5) is converted to The equation of motion for ϕ yields, which, at m = 0, turns out to be the anomaly equation for the axial current j µ 5 =ψγ µ γ 5 ψ. The equation of motion for A 1 yields, By requiring the quantity inside the parenthesis to vanish at spatial infinity, we find Using the equation of motion, or equivalently integrating out the gauge field, we get the pure bosonic action, For 4πm > g/( √ πc), at θ = π there exist the degenerate ground states for ϕ satisfying, Given the two solutions ϕ = (−π ± E)/(2 √ π) with some value E, we find the associated total electric fields respectively. For m g, E approaches π and thus the total electric field E approximately takes the values ±g/2.

Appendix B: Finite-size scaling
We are interested in computing the static susceptibility in the Schwinger model, where χ(k, ω; m) is the retarded Green's function of electric fields. In 1+1 dimensional systems, the electric field can be regarded as a topological charge density E = F 01 . In the thermodynamic limit the Green's function is expected to diverge at the quantum critical point. It is, however, not divergent in finite volume with the system size N . Hence, we attempt to identify the critical points, specified by the fermion mass m * , by fitting the regularized Green's functionχ(m; N ).
Around the critical point ∆m := m − m * ∼ 0, the susceptibility behaves as [70], In a finite system of size L = N a, the susceptibility is with the scaling function Q and the correlation length ξ given by The ansatz (B3) is satisfied in the limit N → ∞ by requiring Note the scaling relation can be read off from the relation between the susceptibility and correlation length (Sec.4 in [75]), Combined with (B2) and (B4) we find the scaling relation, The criticality of lattice Schwinger model belongs to the same universality class as 1+1D transverse Ising model, which in turn equivalent to 2D classical Ising model. Hence the critical exponents are, Finally, we employ the scaling function of the form, with fitting parameters a, b, c, which are scale independent. Thus, our fitting function is χ(m; N ) ≈ N γ/ν ac (N |∆m| ν ) b + c 2 = N 7/4 ac (N |∆m|) b + c 2 .

(B10)
Here m c is a N -dependent fitting parameter that we wish to identify.