Generation of ice states through deep reinforcement learning

TL;DR

This paper introduces a deep reinforcement learning approach to generate ice states in the square ice model, learning physical constraints without prior knowledge and outperforming traditional algorithms in sampling tasks.

Contribution

The study develops a novel RL framework that autonomously learns to generate topologically constrained states, applicable to various models with similar constraints.

Findings

RL agent learns ice rule and loop-closing without prior info

Trained policy effectively samples states in MCMC

Outperforms baseline loop algorithm in benchmarks

Abstract

We present a deep reinforcement learning framework where a machine agent is trained to search for a policy to generate a ground state for the square ice model by exploring the physical environment. After training, the agent is capable of proposing a sequence of local moves to achieve the goal. Analysis of the trained policy and the state value function indicates that the ice rule and loop-closing condition are learned without prior knowledge. We test the trained policy as a sampler in the Markov chain Monte Carlo and benchmark against the baseline loop algorithm. This framework can be generalized to other models with topological constraints where generation of constraint-preserving states is difficult.

Tables

Table 1. Table 1: Neural network architectures for the agent.

Name	Layer	Size
Local Network, NN Policy	linear (relu)	$32$
	linear (relu)	$64$
	linear (relu)	$128$
	parallel (relu)	$7$, $1$
Multi-channel Network, CNN Policy	linear (relu)	$32$
	linear (relu)	$64$
	conv (relu)	$32,3\times 3$, stride $2$
	conv (relu)	$16,3\times 3$, stride $2$
	concat	$1088$
	linear (relu)	$128$
	parallel (relu)	$7$, $1$

Table 2. Table 2: Hyperparameters used in training.

Type	Hyperparameter	Value
Training	Solver type	ADAM
	Learning rate	${10}^{-4}$
	${\beta }_1,{\beta }_2$ (hparam for ADAM)	$0.9,0.999$ (default)
	Total steps	${10}^9$
	Update steps	$120$
	Entropy regularization, $\lambda$	$0.0001$
	Exploration	$ϵ-$greedy, $1\%$
	Gradient clip	$40.0$
	Discount factor $\gamma$	$0.99$
System	Number of parameter server	1 (host)
	Number of workers	8
	Number of policy monitor	1
	TensorBoard process	1
	Evaluation duration	60 seconds
Environment	Lattice size	$L=16$
	Observation shape	$O_g:32\times 32\times 1$, $O_l:10$
	Action shape	$7$
	Stepwise reward	$+0.001$
	Target reward	$+1.0$
	Timeout steps	$1024$

Equations

$$p(s_0\to s_T)=p(s_0)\prod _{t=0}^{T-1}\sum _{a_t}{\pi }_{\theta }(a_t|s_t)p(s_{t+1}|s_t,a_t)$$

$$J({\pi }_{\theta })={\mathbb{E}}_{\pi }[R(\tau )]=\sum _{\tau }{\pi }_{\theta }(\tau )R(\tau ),$$

$${\nabla }_{\theta }J({\pi }_{\theta })={\mathbb{E}}_{\pi }[{\nabla }_{\theta }\log {\pi }_{\theta }(\tau )R(\tau )].$$

$$V^{\pi }(s)={\mathbb{E}}_{\pi }\left[\sum _{k=0}{\gamma }^k{r}_{t+k+1}|s_t=s\right],$$

$${\nabla }_{\theta }J({\pi }_{\theta })={\mathbb{E}}_{\theta }\left[\sum _t{\nabla }_{\theta }\log {\pi }_{\theta }(a_t|s_t)A_{\varphi }^{\pi }(s_t,a_t)\right],$$

$$L(\theta,\phi)=-\left[\sum_{t}\Big{(}\log\pi_{\theta}(a_{t},s_{t})A^{\pi}_{\phi}(t)-||V^{\pi}_{\phi}(s_{t})-R_{t}||^{2}\Big{)}+\lambda\pi_{\theta}(a_{t}|s_{t})\log\pi_{\theta}(a_{t}|s_{t})\right].$$

$$\displaystyle r_{T}(s,a)=\left\{\begin{array}[]{ll}r_{g}=1\text{, if the final state is an ice state}\\ r_{f}=0\text{, otherwise}\\ \end{array}\right.$$

$$S(\mathbf{q})$$

Full text

I Introduction

Ice models Lieb (1967) are the simplest models that can describe the statistical properties of proton arrangement in water ice Pauling (1935), and spin configurations of spin ices Bramwell and Gingras (2001). The ground states in the ice model follow a local constraint, called the ice rule. Due to this constraint, performing local changes in a given ice configuration does not create a new ice state. The loop algorithm, widely used in the Monte Carlo sampling of ice models, generates ice states by exploiting the fact that the difference between two ice configurations are in the form of closed loops Rahman and Stillinger (1972); Yanagawa and Nagle (1979); Barkema and Newman (1998). This is an example where an efficient algorithm generating new constraint-satisfying states emergences from clever inspection of sample states. We want to explore how to automate this type of discovery with machine learning techniques and provide a proof-of-principle implementation.

Efficient state generation is important, in particular, in Markov chain Monte Carlo (MCMC) Landau and Binder (2014); Newman and Barkema (1999). To generate statistically sound samples in the simulation, we need to reduce the autocorrelation between samples in the Markov chain. This can be achieved by clever design of update proposals; for example, in the case of two-dimensional ferromagnetic Ising model, cluster updates such as Swendsen-Wang and Wolf algorithms Swendsen and Wang (1987); Wolff (1989) are powerful methods to reduce the effects of critical slowing down near the critical point and allow for precise extraction of physics with larger systems. However, these cluster algorithms are normally designed for specific models and it is hard to generalize.

Recently, a new class of proposals termed self-learning Monte Carlo, which uses configurations generated by a small size simulation to train an effective model, and draws samples from the effective model to generate new configurations, have demonstrated success in improving of the efficiency of both classical and quantum Monte Carlo simulations Liu et al. (2017a); Nagai et al. (2017); Liu et al. (2017b); Xu et al. (2017). However, the amount of training data required for a successful training normally scales with the complexity of the effective model. Another class of proposals uses the generative models based on deep neural networks in order to generate new samples. The restrict Boltzmann machine, a neural network architecture connected to the the real-space renormalization group Mehta and Schwab (2014), has been used to generate Monte Carlo samples to study thermodynamics Torlai and Melko (2016) and to accelerate Monte Carlo simulations Huang and Wang (2017). The generative adversarial network Goodfellow et al. (2014), where two neural networks contesting with each other, is used to generate configurations for two-dimensional Ising model Liu et al. (2017) and for a complex scalar field in two dimensions with a finite chemical potential Zhou et al. (2018). The success behind these methods relies on the fact that the underlying distribution of the configurations is continuous in the thermodynamic limit; that is, all configurations generated by the neural networks are allowed. On the other hand, direct application of these methods to generate configurations for the topologically constrained model can be difficult. Here not all configurations generated by the neural network can satisfy the constraint. Therefore, a new scheme is highly coveted.

In this paper, we explore how to apply a reinforcement learning framework Sutton and Barto (2018) to generate new configurations that satisfy the ice rule constraint. Reinforcement learning has demonstrated remarkable abilities in achieving better than human performances on video games Mnih et al. (2015), and the game of Go Silver et al. (2016, 2017). Recently it has also been applied to quantum physics research such as quantum control Fösel et al. (2018), quantum error correction Sweke et al. (2018); Liu and Poulin (2018); Poulsen Nautrup et al. (2018), quantum experiment design Melnikov et al. (2018). In a reinforcement learning setup, a machine agent interacts with its environment. At each time step, the agent takes an action on the current state of the environment, and then receives a calculated reward from the observation of the consequence of this action. A feasible reinforcement learning algorithm seeks to maximize the total reward through a trial-and-error learning procedure. We will utilize this feature to design an algorithm to adaptively search for a move from an ice state to another.Similar ideas have been explored recently by the policy-guided Monte Carlo (PGMC) proposal Bojesen (2018) to improve MCMC sampling. The policies in Ref. Bojesen (2018) are modeled with simple functions and the best policy is obtained by optimizing parameters in these models. Here we instead model the policy as a deep neural network, and build the physical constraints into the design of the reward function. Therefore, we can take advantage of the recent progress in the reinforcement learning research to train an agent with a simple set of actions to explore the physical environment. The best policy naturally emerges out of the agent’s interaction with the environment, subjected to the reward design. This scheme should serve as a general framework to search for methods to generate topologically constrained states in other physical models.

The rest of the paper is organized as follows: In Sec. II, we introduce the square ice model and the conventional loop algorithm. In Secs. III and IV, we introduce the deep reinforcement learning framework we use. In Sec. V, the training process and agent behavior will be discussed. We use the trained policies as samplers for Markov chain Monte Carlo in Sec. VI. Finally, we conclude in Sec. VII.

II Square Ice Model

Here we consider a two-dimensional ice model on a square lattice Lieb (1967). At each vertex, the arrangement of spins satisfies the ice rule: two spins pointing into and two pointing out of the vertex (Fig. 1(a)). There exist six types of vertices that satisfy this constraint and all other types of vertices are considered defects with higher energy. This is the reason why ice models are sometimes referred as six-vertex models Lieb (1967); Baxter (2016). Ice states, therefore, correspond to spin configurations with zero defects. All ice states have the same energy, and they form an extensively degenerate ground state manifold with a finite residual entropy at zero temperature of $S/N=k_{B}\ln W$ where $W=\left(4/3\right)^{3/2}$ Lieb (1967). They are, however, sparsely populated in the configuration space of all possible spin configurations (Fig. 2). Consider, for example, flipping a random spin from a given ice configuration, and the two neighboring vertices connected by this spin violates the ice rule constraint, creating a pair of topological defects of three-in-one-out and one-in-three-out vertices (Fig. 1(b)) . The configurations that satisfy the ice rules are thus separated by large energy barriers and the conventional Metropolis single-spin flip update scheme will fail. However, it is possible to flip a series of neighboring spins to bring one ice configuration to another. By examining the difference between two ice states, it is clear that they differ only by spins forming closed-loops. This indicates loop-like updates can bypass the energy barrier and reach a new ice configuration.

The simplest algorithm to generate random closed loops is the long loop algorithm Rahman and Stillinger (1972); Yanagawa and Nagle (1979). In each round of the algorithm, one vertex is chosen as a starting point. One of the two out-going spins from this vertex is chosen with equal probability Barkema and Newman (1998). Following this spin to the connecting vertex and choose one of the two out-going spins with equal probability. Repeat the process until the starting vertex is reached and the loop closes. All the spins on the resulting loop are reversed to update the configuration while preserving the ice rule. Since at each vertex out-going spins are chosen with equal probability, there are $2^{m}$ possible paths for a given starting point, where $m$ is the size of the loop. This algorithm gives an example of how one can generate a new ice configuration by going through a series of highly unfavorable configurations in order to achieve the final goal. This philosophy is indeed very similar to that of the reinforcement learning that the best long-term strategy may involve short-term sacrifices. It is natural to explore the reinforcement learning scheme to generate new ice configurations.

III Reinforcement Learning

Reinforcement learning (RL) is the branch of machine learning concerned with making sequence of decisions for long-term profit  Sutton and Barto (2018). The central idea of RL is to consider a machine agent situated in an environment. At each time step, the agent takes an action on the current state of the environment, and then receives a calculated reward from the observation to the outcome of this action. A feasible RL algorithm seeks to maximize the total reward for the agent from an unknown environment through a trial-and-error learning procedure. We model the agent policy $\pi_{\theta}$ with parameters $\theta$ using deep and convolutional neural networks.

The core idea of our work is to parameterize the proposal operator as the agent policy and search for efficient transitions in the configuration space under given physical constraints. In order to automate the search process, we extend the original Markov chain to a Markov decision process (MDP), by indicating how a state $s_{t}$ transitions into a new state $s_{t+1}$ using the action $a_{t}$ with the transition probability $p(s_{t+1}|s_{t},a_{t})$. The parametrized policy $\pi_{\theta}(a|s)$ specifies the probability of the action $a$ the agent will take given the state $s$. A reward function $r_{t}\equiv r(s_{t},a_{t})$ is associated with each state-action pair $\{s_{t},a_{t}\}$. A move in the Markov chain can, therefore, be decomposed into a series of decision-makings in the MDP (Fig. 3),

[TABLE]

The policy $\pi_{\theta}$ can be regarded as a proposal operator in the conventional MCMC.

We briefly sketch the policy gradient scheme used in this work and refer the interested readers to Refs. Sutton and Barto (2018); Neftci and Averbeck (2019) for details.

Given a trajectory of the state-action pairs $\tau=(\{s_{0},a_{0}\},\{s_{2},a_{2}\},\ldots)$ in the MDP, the goal of RL is to maximize the expected future return, i.e., the weighted sum of rewards. The objective function can be defined as a function of the policy,

[TABLE]

where $\mathbb{E}_{\pi}$ denotes the expectation value is evaluated along a trajectory $\tau$ following the policy $\pi$. The total expected return along the trajectory $\tau$ is given as $R(\tau)=\sum_{t=0}^{T-1}\gamma^{t}r_{t}$, where $0\leq\gamma<1$ denotes a discounting factor. The general method for estimating the gradients of an expectation function is using the score function estimator Shapiro (2003), and the model-free policy gradient can be obtained as

[TABLE]

The optimal parameter $\theta$ can be obtained by the gradient ascent $\theta\leftarrow\theta+\alpha\nabla J(\pi_{\theta})$ with the learning rate $\alpha$.

To estimate how good it is for the agent to be in state $s$, following policy $\pi$, we define a state-value function,

[TABLE]

which estimates the expected future rewards starting from state $s$. To estimate $V^{\pi}$, one starts from a random initial state and follow $\pi(a|s)$ to generate a trajectory to accumulate average return achieved from each state. Given enough samples, the average should converge to the optimal $V_{\pi}$. The value function allows us to compare the policies: a policy $\pi$ is better than or equal to a policy $\pi^{\prime}$ if only if $V_{\pi}(s)\geq V_{\pi^{\prime}}(s)$ for all $s$ Sutton and Barto (2018).

We use the Asynchronous Advantage Actor-Critic (A3C) algorithm Mnih et al. (2016) to perform the training. The value function $V^{\pi}_{\phi}(s)$ is also parametrized by a neural network with parameter $\phi$, and the policy gradient is given by

[TABLE]

where $A^{\pi}_{\phi}(s_{t},a_{t})=r_{t}+\gamma V^{\pi}_{\phi}(s_{t+1})-V^{\pi}_{\phi}(s_{t})$ is an estimate of the advantage function Mnih et al. (2016). The parameterized value function $V^{\pi}_{\phi}(s)$ is trained to predict the expected future return that reduces the variance of policy gradient and stabilizes the training process. In practice, the training is done using the loss function,

[TABLE]

The first term describes the policy loss, the second term gives the estimation of the value function, and the last term provides the entropy regularization. In our experience that the policy entropy term not only promotes exploration but provides a sensible indicator in the training process. The asynchronous algorithm launches several independent agents sampling trajectories from distinct environments. A central network holds shared parameters and is updated asynchronously by many agents. In other words, each agent explores the state space simultaneously from different initial states. This approach prevents the machine from overfitting to some specific configuration and generalizes the policy search. Eight independent agents are used in our program to update eight independent ice states simultaneously.

IV Agent-Environment Interface

We implement an interactive interface for the interaction between the machine agent and the physical environment compatible with OpenAI Gym Brockman et al. (2016). Two kinds of observations of the environment are provided to the agent. The global observation $O_{\rm global}$ monitors the configuration changes $C_{t}=s_{t}-s_{0}$ between the initial ice state $s_{0}$ and the current state $s_{t}$. The local observation $O_{\rm local}$ monitors local information such as the neighboring spins $\sigma_{i}$, the fraction of configuration change $\Delta C$ defined as the ratio of the number of flipped spins and the total number of spins. A state in the MDP is a composite of local and global observations, $s=[O_{\rm local},O_{\rm global}]$.

In deep reinforcement learning, the policy and value functions are both approximated by neural networks. In practice, the value function usually shares weights with the policy approximator Vinyals et al. (2017); therefore, the same network is used for making action decision and predicting the future return. We model the agent with two types of networks: (a) Local network, where only local observation is used as input information; (b) Multi-channel network, where local and global observations are combined. Figure 4 shows the architecture of the multi-channel neural networks for the agent Vinyals et al. (2017). The local information is fed into two layers of feedforward networks while the global observation is extracted through two layers of convolutional neural networks before they are concatenated. The specification of the neural network architecture and hyperparameters for training can be found in App. B. In the following, we call the policy obtained from the local network as NN policy, and multi-channel network as CNN policy.

The actions that the agent located at spin $\sigma$ can execute on the environment contains two types of operations: Select and flip a neighboring spin $\sigma_{i},i\in\{0,..,5\}$, $a_{i}$, and propose an update, $a_{\text{update}}$ (Fig. 6). At each step, the agent is located at a spin and execute one action from the action space according to the policy $\pi_{\theta}(a|s)=[p(a_{0}),\ldots,p(a_{5}),p(a_{\textrm{update}})]$, the probability distribution of taking an action $a$ for a given state $s$. We note that all six neighbors are available in our RL framework, while in the conventional loop algorithm only two are allowed. For example, the policy for the loop algorithm corresponding to the local environment shown in Fig. 6 is $\pi_{\text{loop}}=[0,0,0,0,1/2,1/2,0]$ since only the two outward pointing spins can be chosen with equal probability. Therefore, we design the actions in a less restrictive way and allow the agent to learn the effective trajectory by exploring the physical model. In addition, the action $a_{\text{update}}$ provides flexibility for the agent to terminate at any step. As we will show in the following, after sufficient exploration of the physical environment, the agent tends to perform $a_{\text{update}}$ only when a loop is formed.

The objective of RL is maximizing the cumulative rewards. In playing Atari games, the game score is an obvious choice for reward. In our application, however, searching for update patterns hidden in the environment does not translate directly into an optimization problem. Therefore, the reward function serves as the guiding principal to inform the agent whether a new ice state is generated. The restricted environment itself would be the main reason that versatile policies emerge during the training process Heess et al. (2017). We design the sparse reward function as a step function given at the end of each episode,

[TABLE]

Once the action $a_{\text{update}}$ is executed, the episode will be terminated and the agent will receive the reward $r_{T}$. This encourages the agent to perform $a_{\text{update}}$ only when the environment is in a defect-vacuum state and thus generate a closed loop in each episode. Since the closed-loop condition is relatively rare, at most of the steps the agent receives no reward signals. This kind of approximately silent environment usually causes inefficient sampling and vanishing gradient in the training process. To avoid this problem, we also define a stepwise reward $r_{s}(s,a)$ to encourage the agent to explore the environment and attempt to generate larger loops. The stepwise rewards are assigned at each step regardless of the action taken by the agent and chosen to be relative small such that the agent recognizes $r_{g}$ as the ultimate goal Ng et al. (1999). In practice, we keep the ratio of the stepwise reward to the goal reward as $r_{s}/r_{g}\sim\mathcal{O}(N^{-1})$, where $N$ is the number of sites of the system.

V Training Process and Agent Behavior

At the beginning of the training process, the agent lacks any prior knowledge about the physical environment and performs random moves to collect data. After a few trials, the agent quickly gains the ability to generate a $2\times 2$ square loop, the smallest closed loop providing a non-zero goal reward. In the absence of stepwise reward, the agent would exploit this policy by proposing square loops frequently. However, the stepwise reward encourages the agent to explore the possibilities beyond the square loop. The distribution of loop length as a function of training time is shown in Fig. 7. After about $3\times 10^{5}$ training episodes, the agent becomes proficient at performing larger random loops and executes the update action correctly. While the maximum loop length shows a clear jump at this point, the average loop size grows steadily. The acceptance ratio of the proposed moves also grows steadily during training, and approaches 100% at later stage of the training.

In order to generate a new ice configuration, the agent needs to learn to annihilate the pair of defects by closing the loop. The well-trained agent can be regarded as a policy-guided machine in which the policy $\pi(a|s)$ does the decision-making of local spin flipping and the state-value function $V(s)$ serves as the loop-closing detector. In Fig. 8, we present the decision-making process by the agent near the end of one episode. In most of the steps, the policy distribution shows two possible candidate spins with almost equal probability. This behavior is similar to that of the conventional loop algorithm where one of the two out-going spins are selected with equal probabilities Barkema and Newman (1998). However, there exist time steps that the agent takes deterministic action to move in certain direction. This indicates the agent is exploiting knowledge acquired during the training. Once the loop is closed, the agent decides to perform an update action ($a_{\text{update}}$) and terminate the episode with certainty. This behavior is consistent with the loop closing condition in the loop algorithm. We can also use the state-value function $V(s)$ as an indicator of the agent behavior since it implies the cumulative future reward. In Fig. 9 we show the evolution of $V(s)$ in the process of a single loop formation. When the episode starts, the system is close to an ice state and the value function gives high expectation value. After a few steps of executing the local policy, the value function drops to a local minimum as ends of the segment move apart and going toward the opposite directions (inset A in Fig. 9). During the exploration stage, it is possible that two ends move closer, and the value function reaches a local maximum (inset B in Fig. 9). At the final step, the value function jumps suddenly, anticipating the goal reward to be obtained by performing the update action (inset C in Fig. 9). This indicates that the state-value function plays the role of loop-closing detector and recognizes the global pattern of the environment.

VI Sampling using trained policies

As mention in Sec. IV, two types of policies are trained by supplying different sets of observations to the agent. The NN policy is trained with only the local observation $O_{\textrm{local}}$, and the CNN policy is trained with both the local and global observations $O_{\textrm{local}}$ and $O_{\textrm{global}}$. We will use these trained policies as samplers in the MCMC for the square ice model. With the capacity of deep learning models, it is possible the network simply memorizes the configurations it generated in the learning process. In order to check this, we start the trained agents from the same ice configuration at the same initial position to perform one update. We expect that the agent’s behavior should become independent of the initial starting position and explore the space without special preferences. The baseline loop algorithm shows the probability of each site being visited rapidly smears out (Fig. 10(a)) with the number of steps. We observe, on the other hand, strong memory effects in the NN policy (Fig. 10(b)) such that the agent tends to linger around the initial site; while the CNN policy behaves similar to the baseline loop algorithm (Fig. 10(c)).

We now use the trained policies to propose samples in the Markov chain Monte Carlo to compute physical observables. In particular, we are interested in correlation functions such as the magnetic structure factor defined as

[TABLE]

where the Ising variables $s_{\mathbf{r}}$ on a checkerboard lattice are obtained by an equivalence mapping from the vertex spins on a square lattice (Fig. 11), and $\mathbf{q}$ is the momentum. In order to reproduce the correct correlation, the configurations should be sampled according to the underlying distribution. All ice configurations in the square ice model are equally probable; that is, all ice configurations should be sampled equally for the sampling to be ergodic. Again, we use the structure factor generated by the loop algorithm as our baseline. Fig. 12 shows the spin structure factors using samples generated by different algorithms. Both the NN (Fig. 12(b)) and CNN policies (Fig. 12(c)) can capture the diffuse scattering present in the model; however, it is clear that correlation is biased for the trained policies. We plot in Fig. 12(d) the structure factor $S(\mathbf{q})$ along the high symmetry direction in the first Brillouin zone. We find that some weight of the diffusing scattering is missing in the trained policies, and spectral weight is shifted toward the $\mathbf{Q}=(\pi/2,\pi/2)$ point, indicating bias toward the Néel ordered states for the trained policies. Surprisingly the CNN policy is biased more toward the Néel states than the NN policy, although the former has less memory effects due to the global view. In the Néel order, spins arrange in a regular vortex pattern (Fig. 11). This pattern requires a global view to establish; thus, the CNN policy tends to bias more toward the Néel states.

Clearly, the trained policies are not sampling all the ice states equally, and the ergodicity is broken. This is related to the fact that the machine policy contains actions that makes deterministic moves. These actions indicate that agent somehow memorizes special spin configurations during the training process. Therefore, directly using the policies as samplers introduces bias in the sampling. This is the consequence of not enforcing the detailed balance condition during the training, but may be corrected if one obtains an estimate of the density of the bias in the transition probability. We note this should not be considered a failure of the current framework since the probability distribution of ice states is multimodal with equal weights, which posts a significant challenge to properly sample. On the other hand, by adjusting the stepwise reward, it is possible to bias the agent to create larger loops. Therefore, this can be used to as an overrelaxation step to decorrelate samples in the Markov chain, when combined with the conventional loop algorithm.

VII Conclusion

In this work, we develop a framework using deep reinforcement learning to generate topologically constrained ice states. We successfully train an agent to propose stepwise actions that when combined together can transit from one ice state to another. By rewarding the agent when a correct ice state is generated, several physical insights emerge out of the training. Although the ice rule is not given explicitly, the agent learns to move in a path without creating additional defects. The agent also learns to distinguish open and closed loops, and perform the update action when the loop closes. The physical constraints are built into the reward function, not the policy itself. By using deep reinforcement learning, we show that the machine can actually learn the global loop pattern and propose updates without prior knowledge of the ice rule. Therefore, it is possible to extend this framework to other physical models. For example, quantum Monte Carlo simulations of quantum spin ice models Henry and Roscilde (2014); Wu et al. (2019); Wang et al. (2018), toric code and related gauge models Trebst et al. (2007); Tupitsyn et al. (2010); Kamiya et al. (2015) have been difficult due to the physical constraints present in the model. The generality of the RL framework presented here can explore an enlarged state space and potentially discover new sampling schemes through the automatic exploration of the machine agent on a constrained model. These update proposals can be used as candidate policies in the PGMC Bojesen (2018) where detailed balance can be satisfied. It would also be interesting to model the policy using neural networks with normalizing flows Dinh et al. (2014); Jimenez Rezende and Mohamed (2015); Dinh et al. (2016); Papamakarios et al. (2017); Li and Wang (2018). Combining the exploration capability of RL with the unbiased control of the training and inference of the reversible policy model can potentially lead to the discovery of efficient MCMC algorithms without detailed balance Suwa and Todo (2010) for constrained models. Finally, we note that the training of the model takes roughly three days on a Xeon workstation, and each inference step of the CNN policy takes about $0.75$ ms. Both the conventional loop algorithm and the CNN policy take about five minutes to generate 1000 loops on a $L=32$ system.

VIII Acknowledgments

Y.-J.K. thanks T. Bojesen for useful discussion and critical reading of the manuscript. This work was supported in part by the Ministry of Science and Technology (MOST) of Taiwan under Grants No. 105-2112-M-002-023-MY3, and 107-2112-M-002 -016 -MY3 , 108-2918-I-002 -032, and by the National Science Foundation under Grant No. NSF PHY-1748958. We are also grateful to the National Center for High-performance Computing for computer time and facilities.

Appendix A Sliding horizon scheme

We train the agent policy in a $L_{t}\times L_{t}$ square lattice with the periodic boundary condition. In general, the trained network can only process fixed input dimension assigned before training. In order to use the agent in a larger system, we exploit the translational invariance of the policy and use the sliding horizon embedding that dynamically crops the input image and moves the horizon together with the agent. The observation information within the $L_{t}\times L_{t}$ window is provided to the agent. When the agent reaches the boundary of the window, the window is recentered around the agent as the agent traversing the $L\times L$ lattice. This scheme provides the agent the view similar to the training stage. It also allows us to access system sizes larger than the training size as presented in the main text.

Appendix B Network architecture and hyperparameters

Here we present details of the specifications of the neural network architecture and the hyperparameters used in this paper. Table 1 lists the specification of the local and multi-channel networks for the agent. Table 2 lists all of hyperparameters, system settings and environment configurations for future experiment reproductions.