Tensor network approach to real-time path integral

TL;DR

This paper introduces a tensor network method to represent the real-time path integral for a scalar field in 1+1 dimensions, enabling efficient numerical analysis of quantum field dynamics.

Contribution

The paper develops a tensor network framework for the real-time path integral in scalar field theory, providing a new computational approach for quantum field simulations.

Findings

Tensor network representation matches exact results in non-interacting case

Numerical verification confirms the method's accuracy

Provides a foundation for studying interacting theories

Abstract

We present a tensor network representation of the path integral for the one-component real scalar field theory in 1+1 dimensional Minkowski space-time. It is numerically verified by comparing with the exact result in the non-interacting case.

Figures

Equations

$$S[\varphi ]$$

$$V(\varphi )=\frac{1}{2}(m_0^2-i\epsilon ){\varphi }_x^2+\frac{\lambda }{4!}{\varphi }_x^4,$$

$$\mathcal{Z}$$

$$H^{(0)}(\varphi ,{\varphi }^{\prime })$$

$$H^{(1)}(\varphi ,{\varphi }^{\prime })$$

$${\int }_{-\infty }^{\infty }d\varphi d{\varphi }^{\prime }|H^{(\mu )}(\varphi ,{\varphi }^{\prime }){|}^2<\infty \text{mbox}for\mu =0,1.$$

$$H^{(\mu )}(\varphi ,{\varphi }^{\prime })=\sum _{k=0}^{\infty }{\Psi }_k^{(\mu )}(\varphi ){\sigma }_k{\Phi }_k^{(\mu )\ast }({\varphi }^{\prime }),$$

$${\int }_{-\infty }^{\infty }d\varphi {\Psi }_m^{\ast }(\varphi ){\Psi }_n(\varphi )$$

$$\sum _{k=0}^{\infty }{\Psi }_k(\varphi ){\Psi }_k^{\ast }({\varphi }^{\prime })$$

$$H^{(\mu )}(\varphi ,{\varphi }^{\prime })=\sqrt{2\pi }\sum _{m,n=0}^{\infty }{\psi }_m(\varphi )X_{mn}^{(\mu )}{\psi }_n({\varphi }^{\prime })$$

$${\psi }_n(x)=\frac{1}{\sqrt{{\pi }^{1/2}n!2^n}}{H}_n(x)e^{-x^2/2},H_n(x):\text{mbox}Hermitepolynomials,$$

$$X_{mn}^{(\mu )}=\frac{1}{\sqrt{2\pi }}{\int }_{-\infty }^{\infty }d\varphi d{\varphi }^{\prime }{\psi }_m(\varphi )H^{(\mu )}(\varphi ,{\varphi }^{\prime }){\psi }_n({\varphi }^{\prime }).$$

$$e^{\mp i\varphi {\varphi }^{\prime }}=\sqrt{2\pi }\sum _{n=0}^{\infty }(\mp i{)}^n{\psi }_n(\varphi ){\psi }_n({\varphi }^{\prime }).$$

$$X_{mn}^{(\mu )}=\sum _{k=0}^{\infty }(-1{)}^{k{\delta }_{\mu 0}}{i}^k{G}_{mk}{G}_{nk}\text{mbox}with{G}_{mn}={\int }_{-\infty }^{\infty }d\varphi {\psi }_m(\varphi ){\psi }_n(\varphi )\exp \left[-\frac{i}{4}V(\varphi )\right].$$

$$X_{mn}^{(\mu )}=\sum _{k=0}^{\infty }{U}_{mk}^{(\mu )}{\sigma }_k(V^{(\mu )†}{)}_{kn},\text{mbox}with{U}^{(\mu )},V^{(\mu )}\text{mbox}:unitarymatrix.$$

$$H^{(\mu )}(\varphi ,{\varphi }^{\prime })$$

$${\Psi }_k^{(\mu )}(\varphi )$$

$$T_{ijkl}$$

$$\delta =\left|\frac{\mathrm{Re}[\ln {\mathcal{Z}}_{\mathrm{exact}}]-\mathrm{Re}[\ln {\mathcal{Z}}_N]}{\mathrm{Re}[\ln {\mathcal{Z}}_{\mathrm{exact}}]}\right|.$$

Full text

1 Introduction

The real-time dynamics is one of the most important topics in modern physics, however, the field theoretical approaches based on the path integral representation suffer from the severe sign problem. So far several approaches have been proposed to overcome the sign problem. For example, the complex Langevin method is used to study $\phi^{4}$ scalar field theory [1] and the SU(2) gauge theory [2] in 3+1 dimensional Minkowski space-time. In the course of studies, however, it turns out that the convergence issue gets harder when the time extent is taken relatively larger. Recently, the generalized Lefschetz thimble algorithm is used in the study of the 1+1 dimensional $\phi^{4}$ theory [5] with the Schwinger-Keldysh setup [3, 4]. This approach, however, also has a limitation on the lattice size and it seems hard to perform larger volume simulations.

Tensor network method can also avoid the sign problem and is another promising approach. The idea of tensor networks was invented in the condensate matter physics and developed with a help of quantum information [6, 7, 8, 9]. Recently, the idea and technique have been used in elementary particle physics. It is known that within a framework of the tensor networks there are two types of approaches: the Hamiltonian approach and the Lagrangian one. The former one uses a tensor network in an ansatz form of a trial wave function for the variational method, and the real-time dynamics is studied in various two dimensional models [10, 11, 12]. The letter approach is nothing but the path integral and we will use it in the following. The aim of the study here is to demonstrate how to use the tensor network technique when evaluating the real-time path integral. As a first trial, we work on 1+1 dimensional real scalar field theory with the Minkowskian metric. In the following, we assume the lattice units $a=1$.

2 Tensor network representation

We consider 1+1 dimensional square lattice system whose coordinate is denoted by $x=(x_{0},x_{1})$ with $x_{0}=0,1,2,...,T-1$ and $x_{1}=0,1,2,...,L-1$ for time and space respectively. In the following $T=L$ is assumed. The boundary condition for both directions is taken to be periodic just for a simplicity, though this is not a physically meaningful setting. The one-component real scalar fields on the lattice are denoted by $\phi_{x}\in\mathbb{R}$.

The lattice action for the 1+1 dimensional real scalar $\phi^{4}$ theory with the Minkowskian metric is given by

[TABLE]

with the potential term

[TABLE]

where $m_{0},\lambda\in\mathbb{R}$ are the bare parameters and the Feynman prescription parameter $\varepsilon$ ($\varepsilon\in\mathbb{R}$, $\varepsilon>0$) is introduced and this plays an important role when forming the tensor network representation as we will see later. In an actual computation, note that we will have to take $\varepsilon\to 0$ limit for a bare physical quantity.

In order to use the tensor network method, first of all, one has to rewrite the path integral in terms of the tensor network representation. For that purpose, let us rewrite it as a product of local factors,

[TABLE]

where the local factors for the time-direction $H^{(0)}$ and for the space $H^{(1)}$ are explicitly given by111 Note that the phase factors $\exp(i\phi_{x}^{2}/2)$ associated with the kinetic terms has canceled in eq.(4) and (5).

[TABLE]

They are considered as a two-variable function. An important property of the local factors is that thanks to the presence of the Feynman prescription parameter they are Hilbert-Schmidt operator

[TABLE]

This fact immediately means that they are compact operator and can be expanded as222 One can show that the singular value is independent of $\mu=0,1$. Moreover, it can be shown that it is independent of the bare parameters $m_{0}$ and $\lambda$ since they are included in the phase factor in eq.(4) and (5). Therefore the singular value is a function of only $\varepsilon$, while the bare parameter dependence are captured in the basis functions.

[TABLE]

where $\Psi$ and $\Phi$ are orthonormal basis function satisfying

[TABLE]

and $\sigma_{k}$ represents singular value and is non-negative [13, 14].

As a next step let us see how to obtain the basis functions and the singular values. First, we expand $H^{(\mu)}$

[TABLE]

in terms of the Hermite function $\psi_{m}$ which is the eigen-function of the harmonic oscillator

[TABLE]

and satisfies the orthonormal and the completeness relations as in eq.(8) and (9). Basically $H^{(\mu)}$ are an oscillating function thus it is reasonable to use $\psi_{m}$ as a basis since the latter is also such a function. By using the orthonormal property of $\psi_{m}$, the coefficient matrix $X_{mn}^{(\mu)}$, which is a complex valued symmetric matrix, is obtained as

[TABLE]

An evaluation of the double integral for such an oscillating function is rather technically demanding, thus here we use a trick in the following. In order to separate the integration into two parts, we use a formula333 This may be derived from $e^{ixp}=\sqrt{2\pi}\langle x|p\rangle=\sqrt{2\pi}\sum_{n=0}^{\infty}\langle x|n\rangle\langle n|p\rangle$ and the Fourier transformation of the Hermite function.

[TABLE]

By using the formula together with eq.(4) and (5), $X^{(\mu)}_{mn}$ is given by

[TABLE]

The resulting single integration for $G_{mn}$, which is again a complex valued symmetric matrix, can be exactly carried out for the free case and then $G_{mn}$ is given by the hypergeometric function, while for the interacting case it may be estimated by some numerical integration scheme with a suitable deformation of integration path. After obtaining $X_{mn}^{(\mu)}$ numerically, we apply SVD to it,

[TABLE]

By inserting the above equation into eq.(10), $H^{(\mu)}$ are given as

[TABLE]

where in the last step, we have defined new basis functions

[TABLE]

In this way, the basis functions and the singular values are obtained. In an actual calculation, one cannot deal with the infinite range of index for the matrix, thus one has to truncate it and should monitor the truncation errors of some physical quantities. Note that the singular values of $H^{(\mu)}$ are the same as those of $X^{(\mu)}_{mn}$.

Final step is to obtain the tensor. By collecting all building blocks, the tensor is formed as

[TABLE]

Note that the numbering of the basis function turns out to be an index of tensor $i,j,k$ and $l$. Basically, the range of the indices is infinity but to put it on a computer one has to truncate it. The truncation order is denoted as $N$ and the index range is restricted as $0\leq i,j,k,l\leq N$ in the following. The integration in eq.(18) can be carried out exactly and it is given by the hypergeometric function again.

3 Numerical results

The resulting singular values $\sigma_{k}$ are shown in Fig. 1 (left). For large $\epsilon$ where a damping in the integrand gets stronger, the hierarchy is sharper as expected, while for small $\epsilon$, the hierarchy is not clear and the information compression is not well realized.

To see the hierarchy structure of the tensor elements, $T_{iiii}$ are shown in Fig. 1 (right). Basically the size of elements tends to be exponentially small for larger $i$. This behavior allows us to truncate the tensor index.

In order to check whether the initial tensor in eq.(18) is correctly made, we compute the path integral on $2\times 2$ lattice with the periodic boundary condition. We show the real part of $\ln{\cal Z}/4$ in Fig. 2 (left). Here note that we did not use a coarse-graining, thus this is purely a check of the initial tensor itself. The exact result and the numerical results with some truncation orders are consistent with each other at this scale. Right panel of Fig. 2 shows the relative deviation from the exact results

[TABLE]

By increasing $N$, the deviation $\delta$ tends to be smaller. Note that for smaller $\varepsilon$, the deviation is larger since the hierarchy of the singular value is not so good as seen before.

The path integral at larger volumes can be obtained by using the tensor renormalization group algorithm [15]. The results in Fig. 3 show that the coarse-graining step introduces visible errors. This part should be improved by using some sophisticated algorithms in future .

4 Summary and outlook

We have derived the tensor network representation for the scalar field theory with the Minkowskian metric and numerically checked the validity of the formulation for the free case. One of the important point in the procedure is to expand the local factor $H^{(\mu)}$ in terms of the orthonormal basis functions. The another key point is to introduce the Feynman prescription $\epsilon$ that not only provides a damping factor in the the path integral but also dictates the hierarchy of the singular values of $H^{(\mu)}$. Note that although there is no sign problem in the formulation here, another problem emerges as an ”information incompressibility problem” which means that the hierarchy of the singular value tends to be worse for smaller $\epsilon$. Furthermore we have to eliminate the regulator ($\epsilon\to 0$) in the end of the calculation. This is what we have to pay a price in our formulation.

As future outlooks, there are many things to do. For example, when the $\epsilon$ is small, the initial tensor has large truncation errors thus one may have to improve the initial tensor itself. Secondly, when making the initial tensor for the interacting case, one needs an efficient numerical method to evaluate the oscillating integral, say the steepest descent method and so on. Thirdly, instead of the Feynman prescription, one may use the tilted time axis as a regulator. Finally, one may straightforwardly extend to the Euclidean space as well as the Schwinger-Keldysh formulation in a similar way. Once such a setup is formulated, the transportation coefficients will be accessible.

This work is supported in part by JSPS KAKENHI Grant Numbers JP17K05411 and MEXT as “Exploratory Challenge on Post-K computer” (Frontiers of Basic Science: Challenging the Limits).