Conditional Normalizing flow for Monte Carlo sampling in lattice scalar field theory

In lattice field theory, Monte Carlo Simulation techniques are used to sample lattice configurations based on a distribution defined by the action of the lattice theory. The parameter value at which we generate the lattice samples determines the cost of the simulation. The non-critical region of the lattice theory has low simulation costs for algorithms like Hybrid Monte Carlo (HMC))[1]. However, as we attempt to sample uncorrelated lattice configurations from the critical region, the simulation cost increases rapidly. In the critical region, the integrated autocorrelation time, which gives the measure of correlation, increases rapidly and diverges at the critical point. For a finite-size lattice critical point corresponds to the peak point of the autocorrelation curve. This problem is known as the critical slowing down[2, 3]. Many efforts have been made to lessen the impact of critical slowing in statistical systems and lattice QFT [4, 5, 6]. But it always remains a challenging task to simulate near the critical point of a lattice QFT by overcoming the critical slowing down.

These days, ML-based solutions to this problem are becoming popular. Various ML algorithms have been applied for statistical physics and condensed matter problems[7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19] . Some generative learning algorithms[21, 22, 23, 24, 25, 26, 27] have recently been developed to avoid the difficulty in lattice field theory. Conditional GAN was used in a 2D lattice Gross Neveue model [20] and shown to be effective in mitigating the influence of critical slowdown. However, explicit probability density estimation is not accessible in GAN; thus, we cannot guarantee that the model distribution is identical to the true lattice distribution. Flow based generative leanings are found successful in avoiding the problem of critical slowing down in scalar field theory[21], fermionic system[22], U(1) gauge theory[24] and Schwinger models[26]. In the flow-based approach [21, 24] an NF model is trained at a single value of the action parameter with reverse KL divergence. Finally, the trained model can generate lattice samples at the same parameter value. This model is initialized from scratch for a parameter value and does not use any lattice samples. Since HMC simulation in the non-critical region is not affected by the critical slowing down problem, we can use samples from that region to train a generative model. We present a method for sampling lattice configurations near the critical regions using Conditional Normalizing Flow (C-NF) to reduce the problem of critical slowing down. This method involves training a C-NF model in a non-critical region to produce samples for several parameter values in the critical region. Our goal is to generate lattice configurations from the distribution $p(\phi|\lambda_{crit})=\frac{1}{Z}e^{-S(\phi|\lambda_{crit})}$, where $\phi$ denotes the lattice field, $\lambda_{crit}$ denotes the action parameter close to critical point and Z is the partition function. We train a C-NF model $\tilde{p}({\phi|\lambda})$ with HMC samples from $p(\phi|\lambda)$ for various non-critical $\lambda$ values. We train the C-NF model to be a generalized model over $\lambda$ parameters. The model is then interpolated or extrapolated to the $\lambda$ values in the critical region to generate lattice configurations. However, the interpolated or extrapolated model may not directly provide samples from the true distribution. But the exactness can be guaranteed by using the Metropolis-Hastings(MH) algorithm at the end. So, after training, we use the interpolated/extrapolated model $\tilde{p}(\phi|{\lambda_{crit}})$ at critical region as proposal for constructing a Markov Chain via an independent MH algorithm[21]. This method is useful when the probability distribution is known up to a normalizing factor.

The primary contributions of this study are as follows:

1. 1.

   Using Conditional Normalizing Flows, we present a new method for sampling lattice configurations near the critical regions. This method eliminates the critical slowdown problem.

2. 2.

   The model has the ability to learn about the lattice system across multiple $\lambda$ values and use this knowledge to generate sample at any given $\lambda$ values. As a result, our model can generate samples at multiple $\lambda$ values values in the critical region, which is not possible for the existing flow-based methods in lattice thoery.

3. 3.

   We also demonstrate that Conditional Normalizing flow can be used to do both interpolation and extrapolation for lattice $\phi^{4}$ theory. The extrapolation demonstrates the possibility of using our approach for sampling lattice gauge theory.

In $2d$ euclidean space, the action for $\phi^{4}$ theory can be written as:

$\displaystyle S(\phi,\lambda,m)=\int dx^{2}[(\partial_{\mu}\phi(x))^{2}+m^{2}\phi(x)^{2}+\lambda\phi(x)^{4}]$

(1)

where $\lambda$ and $m$ are the two parameters of the theory.

On lattice the action become:

	$\displaystyle S(\phi,\lambda,m)=$	$\displaystyle\sum\limits_{x}\Big{[}\sum\limits_{\mu={1,2}}\Big{[}2\phi(x)^{2}+\phi(x)\phi(x+\hat{\mu})$
		$\displaystyle-\phi(x)\phi(x-\hat{\mu})\Big{]}+m^{2}\phi(x)^{2}+\lambda\phi(x)^{4}\Big{]}$		(2)

where $x$ is a $2d$ discrete vector and $\hat{\mu}$ represents two possible directions on the lattice. $\phi(x)$ is defined on each lattice site, taking only real values.

We choose this specific form of action for HMC simulation because it is suited for creating datasets for training the C-NF model. Using this form, we do not need to apply any further transformations to the samples during training. This action possesses $\phi(x)\rightarrow\phi(-x)$ symmetry, however it is spontaneously broken at a specific parameter region. We set $m^{2}=-4$ for our numerical experiments and observed spontaneous symmetry breaking in the theory by varying $\lambda$. If we begin in the broken phase of the theory, the order parameter under consideration $\langle\phi^{2}\rangle$ is nonzero which approaches zero at the critical point and remains zero in the symmetric phase, as shown in Figure 1. In HMC simulation we choose a parameter $\lambda$ and produce configurations based on the probability distribution:

	$\displaystyle P(\phi|\lambda)=\frac{1}{Z}e^{-S(\phi|\lambda)}$		(3)
	$\displaystyle\text{where, }Z=\sum\limits_{\phi}e^{-S(\phi|\lambda)}$

Each lattice configuration is a $2d$ matrix with the dimensions $L_{x}\times L_{y}$. For our experiment, we use $L_{x}=L_{y}=8$ and add a periodic boundary condition to the lattice in all directions. More information on lattice $\phi^{4}$ theory can be found in ref. [28]. Some of the observables which we calculate on the lattice ensembles are:

1. 1.

   $\langle\tilde{\phi}^{2}\rangle$: $\tilde{\phi}=\frac{1}{V}\sum_{x}\phi(x)$

2. 2.

   Correlation Function:

   $\displaystyle G_{c}(x)=\frac{1}{V}\sum_{y}[\langle\phi(y)\phi(x+y)\rangle-\langle\phi(y)\rangle\langle\phi(x+y)\rangle]$

   Zero momentum Correlation Function:
   $C(t)=\sum_{x_{1}}G_{c}(x_{1},t)$

3. 3.

   Two Point Susceptibility: $\chi=\sum_{x}G(x)$

Normalizing flows[29] are a generative model for constructing complex distributions by transforming a simple known distribution via a series of invertible and smooth mapping $f:\mathbb{R}^{d}\rightarrow\mathbb{R}^{d}$ with inverse $f^{-1}=g$. If $p_{z}(z)$ is the prior distribution and $p_{t}(x)$ is the complex target distribution, then the model distribution $p_{m}(x;\theta)$ can be written using the change of variable formula as

	$\displaystyle p_{m}(x;\theta)=p_{z}(z)\big{|}det\frac{\partial f_{\theta}(z)^{-1}}{\partial x}\big{|}$		(4)
	$\displaystyle\text{where, }x=f_{\theta}(z)$

Fitting a flow-based model $p_{m}(x;\theta)$ to a target distribution $p_{t}(x)$ can be accomplished by minimising their KL divergence. The most crucial step is to build the flow so that we can calculate $\big{|}det\frac{\partial f_{\theta}(z)^{-1}}{\partial x}\big{|}$. One such method is the affine coupling block, which divides the input $z$ into two halves and applies an affine transformation to produce upper or lower triangular Jacobians. The transformation rules for such a building are as follows[30]:

	$\displaystyle x_{1}=z_{1}\odot exp(s_{1}(z_{2})+t_{1}(z_{2}))$
	$\displaystyle x_{2}=z_{2}\odot exp(s_{2}(x_{1}))+t_{2}(x_{1})$		(5)

where, $\odot$ represent element-wise product of two vectors.

The inverse of this coupling layer is simply computed as:

	$\displaystyle z_{2}=(x_{2}-t_{2}(x_{1}))\odot exp(-s_{2}(x_{1}))$
	$\displaystyle z_{1}=(x_{1}-t_{1}(z_{2}))\odot exp(-s_{1}(z_{2}))$		(6)

Because inverting an affine coupling layer does not require the inverse of $s_{1},s_{2},t_{1},t_{2}$, they can be any non-linear complex function and can thus be represented by neural networks. Introducing a conditioning parameter in NF is not as simple as it is in GAN. However, since $s$ and $t$ are only evaluated in the forward direction, we may concatenate the conditioning parameter $c$ with the input $x$ to the coupling layer in order to invert the model[31], as illustrated in Figure 2. As a result, the affine coupling transformation rules become:

	$\displaystyle x_{1}=z_{1}\odot exp(s_{1}(z_{2},c)+t_{1}(z_{2},c))$
	$\displaystyle x_{2}=z_{2}\odot exp(s_{2}(x_{1},c))+t_{2}(x_{1},c))$		(7)

And its inverse become

	$\displaystyle z_{2}=(x_{2}-t_{2}(x_{1},c))\odot exp(-s_{2}(x_{1},c))$
	$\displaystyle z_{1}=(x_{1}-t_{1}(z_{2},c))\odot exp(-s_{1}(z_{2},c))$		(8)

Let us designate the C-NF model as $f(x;c,\theta)$ and its inverse as $g(z;c,\theta)$. The invertibility for any fixed condition $c$ is given by

$\displaystyle f^{-1}(.;c,\theta)=g(.;,c,\theta)$

(9)

The change of variable formula become

$\displaystyle p_{m}(x;c,\theta)=p_{z}(f(x;c,\theta))|detJ_{f}(z)|^{-1}$

(10)

And the loss function is the KL divergence between the model distribution $p_{m}(x;c,\theta)$ and target distribution $p_{t}(x;c)$ :

$\displaystyle\mathcal{L}=D_{KL}[{p_{t}(x;c)||p_{m}(x;c,\theta)}]$

(11)

We can find the maximum likelihood network parameter $\theta_{ML}$ using this loss. Then, for a fixed $c$, we can execute conditional generation of $x$ by sampling $z$ from $p_{z}(z)$ and employing the inverted network $g(z;c,\theta_{ML})$.

This section discusses dataset preparation and the model architecture utilized in training the C-NF model. In addition, the training details of the C-NF model and the sampling process in the critical region are discussed. We train the C-NF model with different datasets and model architecture for interpolation and extrapolation.

The sampling algorithm for the baseline approach (HMC) and the suggested method is vastly different; thus, a direct cost comparison is opaque. Nonetheless, we separate the simulation cost for the proposed technique into two components: Training time and sample generation time. On a Colab Tesla P100 GPU, the training time for the C-NF model for interpolation or extrapolation is roughly 5-6 hours. However, sample generation is very fast for the C-NF model. Generating one Markov Chain of $10^{5}$ configuration takes 5-7 minutes with a 25-40$\%$ acceptance rate. Due to the short generation time, a low acceptance rate is acceptable until the autocorrelation time increases.

Once the C-NF model has been trained, it can be employed repeatedly to generate configurations for a wide range of $\lambda$ values. From a single training of the C-NF(interpolated) model, we have generated configurations for 13 $\lambda$’s in the critical region. So, our approach outperforms HMC for sampling at multiple $\lambda$ values in the critical region.

The critical slowing down problem prevents generating a large ensemble in the critical region of a lattice theory. In order to resolve this, we employ a Conditional normalizing flow trained on HMC samples with low autocorrelation and generate samples in the critical region. In order to learn a general distribution over parameter $\lambda$, we train the C-NF model away from the critical point of lattice $\phi^{4}$ theory. This model is interpolated in the critical region and serves as a proposal for the MH algorithm to generate a Markov chain. The degree to which the extrapolated or interpolated model resembles the actual distribution determines the acceptance rate in the critical region. In order to achieve a high acceptance rate and prevent the development of autocorrelation, the C-NF model must be trained adequately. The C-NF model generates uncorrelated samples, and we trained well enough to get 25-45$\%$ acceptance rate. With this much acceptance rate, we found no correlation between configuration in the Markov chain. As a result, our method significantly mitigates the critical slowing down problem. Aside from that, our method can be highly efficient when we need interpolation/extrapolation to numerous $\lambda$ values in the critical region. Since lattice gauge theory requires extrapolation to the critical region, we have likewise extrapolated $\phi^{4}$ lattice theory to the critical region. We observe high agreement between observables estimated using the suggested technique and HMC simulation for both interpolation and extrapolation.