Reinforcement Learning in Different Phases of Quantum Control

TL;DR

This paper demonstrates that reinforcement learning can effectively find high-fidelity quantum control protocols in complex many-body systems, performing comparably to optimal control methods and revealing a phase transition in protocol space.

Contribution

The study introduces RL techniques for quantum state control, showing their ability to discover near-optimal protocols in complex quantum systems and uncovering a phase transition in control landscape.

Findings

RL achieves high-fidelity control comparable to optimal methods

Quantum control exhibits a spin-glass-like phase transition

RL can find near-optimal protocols even in glassy phases

Abstract

The ability to prepare a physical system in a desired quantum state is central to many areas of physics such as nuclear magnetic resonance, cold atoms, and quantum computing. Yet, preparing states quickly and with high fidelity remains a formidable challenge. In this work we implement cutting-edge Reinforcement Learning (RL) techniques and show that their performance is comparable to optimal control methods in the task of finding short, high-fidelity driving protocol from an initial to a target state in non-integrable many-body quantum systems of interacting qubits. RL methods learn about the underlying physical system solely through a single scalar reward (the fidelity of the resulting state) calculated from numerical simulations of the physical system. We further show that quantum state manipulation, viewed as an optimization problem, exhibits a spin-glass-like phase transition in the space of protocols as a function of the protocol duration. Our RL-aided approach helps identify variational protocols with nearly optimal fidelity, even in the glassy phase, where optimal state manipulation is exponentially hard. This study highlights the potential usefulness of RL for applications in out-of-equilibrium quantum physics.

Equations

$$\mathcal{S}=\{s=(t,h_x(t))\},\mathcal{A}=\{a=\delta h_x\},\mathcal{R}=\{r\in [0,1]\}.$$

$$r(t)=\{\begin{array}{ll}\text{hfill}0,\text{hfill}& \text{ if t<T}\\ \text{hfill}{F}_h(T)=|⟨{\psi }_{\ast }|\psi (T)⟩{|}^2,\text{hfill}& \text{ if t=T}\end{array}$$

$$\mathrm{Environment}$$

$$s_0\to a_0\to r_0⟶s_1\to a_1\to r_1⟶s_2\to \cdots ⟶s_{N_T}.$$

$$Q(s_i,a_i)⟵Q(s_i,a_i)+\alpha \left[r_i+\underset{a}{\max }Q(s_{i+1},a)-Q(s_i,a_i)\right],$$

$$H[h_x(t)]=-S^z-h_x(t)S^x,$$

$$q(T)=\frac{1}{16N_T}\sum _{j=1}^{N_T}\overline{\{h_x(j\delta t)-\overline{h_x}(j\delta t){\}}^2},$$

$$H[h_x(t)]=-\sum _{j=1}^L\left(S_{j+1}^z{S}_j^z+g{S}_j^z+h_x(t)S_j^x\right).$$

$${\mathcal{F}}_h({\tau }^{(1)},T-{\tau }^{(1)})$$

$$H[h_x]$$

$${\mathcal{F}}_h^{2D}({\tau }^{(1)},{\tau }^{(2)},T-{\tau }^{(1)}-{\tau }^{(2)})$$

$$H[h_x]$$

$$h_{\mathrm{CRAB}}(t)$$

$$\mathcal{C}[h(t)]=\mathcal{F}(\{A_i,B_i,{\omega }_i\})+\frac{1}{16T}{\int }_0^T\mathrm{d}t[h(t){]}^2$$

$$F_h(T)=|⟨{\psi }_{\ast }|U(T,0)|{\psi }_i⟩{|}^2,$$

$$\frac{\delta F_h(T)}{\delta h_x(t)}=i⟨{\psi }_{\ast }|U(T,0)|{\psi }_i⟩⟨{\psi }_i|U(0,t)XU(t,T)|{\psi }_{\ast }⟩-i⟨{\psi }_{\ast }|U(T,t)XU(t,0)|{\psi }_i⟩⟨{\psi }_i|U(0,T)|{\psi }_{\ast }⟩$$

$$\frac{\delta F_h(T)}{\delta h_x(t)}=2\mathrm{Im}\left[⟨\varphi (t)|X|{\psi }_i(t)⟩\right],$$

$$H(t)=-S^z-h_x(t)S^x,|{\psi }_i⟩\sim (-1/2-\sqrt{5}/2,1{)}^{\mathsf{T}},|{\psi }_{\ast }⟩\sim (1/2+\sqrt{5}/2,1{)}^{\mathsf{T}},$$

$$h_x^{\mathrm{optimal}}(t)=\underset{h_x(t):|h_x(t)|\le 4}{\mathrm{argmin}}\{I_h(T)\}=\underset{h_x(t):|h_x(t)|\le 4}{\mathrm{argmax}}|⟨{\psi }_{\ast }|{\mathcal{T}}_t{\mathrm{e}}^{-i{\int }_0^T\mathrm{d}tH[h_x(t)]}|{\psi }_i⟩{|}^2,$$

$$q(T)=\frac{1}{16N_T}\overline{\sum _{n=1}^{N_T}\{h_x(n\delta t)-\overline{h_x}(n\delta t){\}}^2},\overline{h_x}(t)=\frac{1}{N_{\mathrm{real}}}\sum _{\alpha =1}^{N_{\mathrm{real}}}{h}_x^{\alpha }(t).$$

Full text

I Introduction

Reliable quantum state manipulation is essential for many areas of physics ranging from nuclear magnetic resonance experiments Vandersypen and Chuang (2005) and cold atomic systems van Frank et al. (2016); Wigley et al. (2016) to trapped ions Islam et al. (2011); Senko et al. (2015); Jurcevic et al. (2014), quantum optics Sayrin et al. (2011), superconducting qubits Barends et al. (2016), nitrogen vacancy centers Zhou et al. (2017), and quantum computing Nielsen and Chuang (2002). However, finding optimal control sequences in such experimental platforms presents a formidable challenge due to our limited theoretical understanding of nonequilibrium quantum systems, and the intrinsic complexity of simulating large quantum many-body systems.

For long protocol durations, adiabatic evolution can be used to robustly reach target quantum states, provided the change in the Hamiltonian is slow compared to the minimum energy gap. Unfortunately, this assumption is often violated in real-life applications. Typical experiments often have stringent constraints on control parameters, such as a maximum magnetic-field strength or a maximal switching frequency. Moreover, decoherence phenomena impose insurmountable time constraints beyond which quantum information is lost irreversibly. For this reason, many experimentally-relevant systems are in practice uncontrollable, i.e. there are no finite-duration protocols, which prepare the desired state with unit fidelity. In fact, in Anderson and many-body localized, or periodically-driven systems, which are naturally away from equilibrium, the adiabatic limit does not even exist Khemani et al. (2015); Weinberg et al. (2017). This has motivated numerous approaches to quantum state control Baksic et al. (2016); Wang et al. (2015); Agundez et al. (2017); Bao et al. (2017); Rotskoff et al. (2017); Leung et al. (2017); Yang et al. (2017); Jarzynski (2013); Kolodrubetz et al. (2017); Sels and Polkovnikov (2017); Glaser et al. (1998); Rabitz et al. (2000); Khaneja et al. (2001); Sklarz and Tannor (2002); Demirplak and Rice (2003); Khaneja et al. (2005a); Berry (2009); Caneva et al. (2011); Zahedinejad et al. (2014, 2015); Theisen et al. (2017); Bukov et al. (2018a); Wu et al. (2018). Despite all advances, at present date surprisingly little is known about how to successfully load a non-integrable interacting quantum system into a desired target state, especially in short times, or even when this is feasible in the first place Caneva et al. (2011); Doria et al. (2011); Caneva et al. (2014); Lloyd and Montangero (2014).

In this paper, we adopt a radically different approach to this problem based on machine learning (ML) Sutton and Barto (2017); Bishop (2007); Mnih et al. (2015); Silver et al. (2016); Dunjko et al. (2016); Mehta et al. (2018); Judson and Rabitz (1992). ML has recently been applied successfully to several problems in equilibrium condensed matter physics Carleo and Troyer (2017); Carrasquilla and Melko (2017), turbulent dynamics Reddy et al. (2016); Colabrese et al. (2017) and experimental design Krenn et al. (2016); Melnikov et al. (2017), and here we demonstrate that Reinforcement Learning (RL) provides deep insights into nonequilibrium quantum dynamics Chen et al. (2014); August and Hernández-Lobato (2018); Fösel et al. (2018); Zhang et al. (2018); Niu et al. (2018); Albarran-Arriagada et al. (2018). Specifically, we use a modified version of the Watkins Q-Learning algorithm Sutton and Barto (2017) to teach a computer agent to find driving protocols which prepare a quantum system in a target state $|\psi_{*}\rangle$ starting from an initial state $|\psi_{i}\rangle$ by controlling a time-dependent field. A far-reaching consequence of our study is the existence of phase transitions in the quantum control landscape of the generic many-body quantum control problem. The glassy nature of the prevalent phase implies that the optimal protocol is exponentially difficult to find. However, as we demonstrate, the optimal solution is unstable to local perturbations. Instead, we discover classes of RL-motivated stable suboptimal protocolsGingrich et al. (2016), the performance of which rival that of the optimal solution. Analyzing these suboptimal protocols, we construct a variational theory, which demonstrates that the behaviour of physical degrees of freedom (d.o.f.) [which scale exponentially with the system size $L$ for ergodic models] in a non-integrable many-body quantum spin chain can be effectively described by only a few variables within the variational theory. We benchmark the RL results using Stochastic Descent (SD), and compare them to optimal control methods such as CRAB Caneva et al. (2011) and (for simplicity) first-order GRAPE Khaneja et al. (2005b) (without its quasi-Newton extensions Machnes et al. (2011); Agundez et al. (2017); De Fouquieres et al. (2011)), see discussion in the Supplemental Material sup .

In stark contrast to most approaches to quantum optimal control, RL is a model-free feedback-based method which could allow for the discovery of controls even when accurate models of the system are unknown, or parameters in the model are uncertain. A potential advantage of RL over traditional derivative-based optimal control approaches is the fine balance between exploitation of already obtained knowledge and exploration in uncharted parts of the control landscape. Below the quantum speed limit Hegerfeldt (2013), exploration becomes vital and offers an alternative to the prevalent paradigm of multi-starting local gradient optimizers Sørensen et al. (2016). Unlike these methods, the RL agent progressively learns to build a model of the optimization landscape in such a way that the protocols it finds are stable to sampling noise. In this regard, RL-based approaches may be particularly well-suited to work with experimental data and do not require explicit knowledge of local gradients of the control landscape Khaneja et al. (2005b); sup . This may offer a considerable advantage in controlling realistic systems where constructing a reliable effective model is infeasible, for example due to disorder or dislocations.

To manipulate the quantum system, our computer agent constructs piecewise-constant protocols of duration $T$ by choosing a drive protocol strength $h_{x}(t)$ at each time $t=j\delta t$, $j=\left\{0,1,\cdots,T/\delta t\right\}$, with $\delta t$ the time-step size. In order to make the agent learn, it is given a reward for every protocol it constructs – the fidelity $F_{h}(T)=|\langle\psi_{\ast}|\psi(T)\rangle|^{2}$ for being in the target state after time $T$ following the protocol $h_{x}(t)$ under unitary Schrödinger evolution. The goal of the agent is to maximize the reward in a series of attempts. Deprived of any knowledge about the underlying physical model, the agent collects information about already tried protocols, based on which it constructs new, improved protocols through a sophisticated biased sampling algorithm. In realistic applications, one does not have access to infinite control fields; for this reason, we restrict to fields $h_{x}(t)\in[-4,4]$, see Fig. 1b. For reasons relating to the simplicity and efficiency of the numerical simulations, throughout this work we further restrict the RL algorithm to the family of bang-bang protocols 111We note that there is no proof that the family of bang-bang protocols contains the optimal protocol for the single-qubit control problem in Sec. III, since the present control problem is of bilinear type Schättler and Ledzewicz (2012) . An additional advantage of focusing on bang-bang protocols is that this allows us to interpret the control phase transitions we find using the language of Statistical Mechanics Day et al. (2018).

II Reinforcement Learning

Reinforcement Learning (RL) is a subfield of Machine Learning (ML) in which a computer agent learns to perform and master a specific task by exerting a series of actions in order to maximize a reward function, as a result of interaction with its environment. Here, we use a modified version of Watkins online, off-policy Q-Learning algorithm with linear function approximation and eligibility traces Sutton and Barto (2017) to teach our RL agent to find protocols of optimal fidelity. Let us we briefly summarize the details of the procedure. For a detailed description of the standard Q-learning algorithm, we refer the reader to Ref. Sutton and Barto (2017).

The fidelity optimization problem is defined as an episodic, undiscounted Reinforcement Learning task. Each episode takes a fixed number of steps $N_{T}=T/\delta t$, where $T$ is the total protocol duration, and $\delta t$ – the physical (protocol) time step. We define the state $\mathcal{S}$, action $\mathcal{A}$ and reward $\mathcal{R}$ spaces, respectively, as

[TABLE]

The state space $\mathcal{S}$ consists of all tuples $(t,h_{x}(t))$ of time $t$ and the corresponding magnetic field $h_{x}(t)$. Notice that with this choice no information about the physical quantum state whatsoever is encoded in the RL state, and hence the RL algorithm is model-free. Thus, the RL agent will be able to learn circumventing the difficulties associated with the theoretical notions in quantum physics. Including time $t$ to the state is not common in Q-Learning, but is required here in order for the agent to be able to estimate how far away it is from the episode’s end, and adjust its actions accordingly. Even though there is only one control field, the space of available protocols grows exponentially with the inverse step size $\delta t^{-1}$.

The action space $\mathcal{A}$ consists of all jumps $\delta h_{x}$ in the protocol $h_{x}(t)$. Thus, protocols are constructed as piecewise-constant functions. We restrict the available actions of the RL agent in every state $s$ such that at all times the field $h_{x}(t)$ is in the interval $[-4,4]$. We verify that RL also works for quasi-continuous protocols with many different steps $\delta h_{x}$ sup . The bang-bang protocols discussed in the next section and the quasi-continuous protocols, used in the Supplemental Material sup , are examples of the family of protocol functions we allow in the simulation.

Last but not least, the reward space $\mathcal{R}$ is the space of all real numbers in the interval $[0,1]$. The rewards for the agent are given only at the end of each episode, according to:

[TABLE]

This reflects the fact that we are not interested in which quantum state the physical system is in during the evolution; all that matters for our purpose is to maximize the final fidelity.

An essential part of setting up the RL problem is to define the environment, with which the agent interacts in order to learn. We choose this to consist of the Schrödinger initial value problem, together with the target state:

[TABLE]

where $H[h_{x}(t)]$ is the Hamiltonian, see Sec. III, whose time dependence is defined through the magnetic filed $h_{x}(t)$ which the agent is constructing during the episode via online Q-Learning updates for specific single-particle and many-body examples.

Let us now briefly illustrate the protocol construction algorithm: for instance, if we start in the initial RL state $s_{0}\!=\!(t\!=\!0,h_{x}\!=\!-4)$, and take the action $a\!=\!\delta h_{x}\!=\!8$, we go to the next RL state $s_{1}=(\delta t,+4)$. As a result of the interaction with the environment, the initial quantum state is evolved forward in time for one time step (from time $t_{0}=0$ to time $t_{1}=\delta t$) with the constant Hamiltonian $H[h_{x}=4]$: $|\psi(\delta t)\rangle=\mathrm{e}^{-iH[h_{x}=4]\delta t}|\psi_{i}\rangle$. After each step we compute the local reward according to Eq. (1), and update the Q-function, even though the instantaneous reward at that step might be zero [the update will still be non-trivial in the later episodes, since information is propagated backwards from the end of the episode, see Eq. (2)]. This procedure is repeated until the end of the episode is reached at $t=T$. In general, one can imagine this partially-observable Markov decision process as a state-action-reward chain

[TABLE]

The above paragraph explains how to choose actions according to some fixed policy $\pi(a|s)$ – the probability of taking the action $a$ from the state $s$. Some RL algorithms, such as Policy Gradient directly optimize the policy. Instead, Watkins Q-Learning offers an alternative which allows to circumvent this. The central object in Q-Learning is the $Q(s,a)$ function which is given by the expected total return $R=\sum_{i=0}^{N_{T}}r_{i}$ at the end of each episode, starting from a fixed state $s$, taking the fixed action $a$, and acting optimally afterwards. Clearly, if we have the optimal $Q$-function $Q_{\ast}$, then the optimal policy is the deterministic policy $\pi_{\ast}(a|s)=1,\ \mathrm{if}\ a=\operatorname*{arg\,max}_{a^{\prime}}Q(s,a^{\prime})$, and $\pi_{\ast}(a|s)=0$ for all other actions.

Hence, in Q-Learning one looks directly for the optimal Q-function. It satisfies the Bellman optimality equation, the solution of which cannot be obtained in a closed form for complex many-body systems 222Bellman’s equation probably admits a closed solution for the single-qubit example from Sec. III.. The underlying reason for this can be traced back to the non-integrability of the dynamical many-body system, as a result of which the solution of the Schrödinger equation cannot be written down as a closed-form expression even for a fixed protocol, and the situation is much more complicated when one starts optimizing over a family of protocols. The usual way of solving the Bellman equation numerically is Temporal Difference learning, which results in the following Q-Learning update rule Sutton and Barto (2017)

[TABLE]

where the learning rate $\alpha\in(0,1)$. Whenever $\alpha\approx 1$, the convergence of the update rule (2) can be slowed down or even precluded, in cases where the Bellman error $\delta_{t}=r_{i}+\max_{a}Q(s_{i+1},a)-Q(s_{i},a_{i})$ becomes significant. On the contrary, $\alpha\approx 0$ corresponds to very slow learning. Thus, the optimal value for the learning rate lies in between, and is determined empirically for the problem under consideration.

To allow for the efficient implementation of piecewise-constant drives, i.e. bang-bang protocols with a large number of bang modes, cf. Ref. sup , we employ a linear function approximation to the $Q$-function, using equally-spaced tilings along the entire range of $h_{x}(t)\in[-4,4]$ Sutton and Barto (2017). The variational parameters of the linear approximator are found iteratively using Gradient Descent. This allows the RL agent to generalize, i.e. gain information about the fidelity of not yet encountered protocols.

We iterate the algorithm for $2\times 10^{4}$ episodes. The exploration-exploitation dilemma Sutton and Barto (2017) requires a fair amount of exploration, in order to ensure that the agent visits large parts of the RL state space which prevents it from getting stuck in a local maximum of reward space from the beginning. Too much exploration, and the agent will not be able to learn. On the other hand, no exploration whatsoever guarantees that the agent will repeat deterministically a given policy, though it will be unclear whether there exists a better, yet unseen one. In the longer run, we cannot preclude the agent from ending up in a local maximum. In such cases, we run the algorithm multiple times starting from a random initial condition, and post-select the outcome. Hence, the RL solution is almost-optimal in the sense that its fidelity is close to the true global optimal fidelity. Unfortunately, the true optimal fidelity for nonintegrable many-body systems is unknown, and it is a definitive feature of glassy landscapes, see Sec. V, that the true optimal is exponentially hard, and therefore also impractical, to find Day et al. (2018).

We also verified that RL does not depend on the initial condition chosen, provided the change is small. For instance, if one chooses different initial and target states which are both paramagnetic, then RL works with marginal drops in fidelity, which depend parametrically on the deviation from the initial and target states. If however, the target is, e.g. paramagnetic and we choose an antiferromagnetic initial state [i.e. the initial and target states are chosen in two different phases of matter], then we observe a drop in the found fidelity.

Due to the extremely large state space, we employ a replay schedule to ensure that our RL algorithm could learn from the high fidelity protocols it encountered. Our replay algorithm alternates between two different ways of training the RL agent which we call training stages: an “exploratory” training stage where the RL agent exploits the current $Q$-function to explore, and a “replay” training stage where we replay the best encountered protocol. This form of replay, to the best of our knowledge, has not been used previously. In the exploratory training stage, which lasts $40$ episodes, the agent takes actions according to a softmax probability distribution based on the instantaneous values of the $Q$-function. In other words, at each time step, the RL agent looks up the instantaneous values $Q(s,:)$ corresponding to all available actions, and computes a probability for each action: $P(a)\sim\exp(\beta_{\mathrm{RL}}Q(s,a))$. This exploration scheme results in random flips in the bangs of the protocol sequence, which is essentially a variation on the instantaneous RL best solution. Fig. 2 shows that some of these variations lead to drastic reduction in fidelity, which we related to the glassy character of the correlated control phase, see Sec. V.

The amount of exploration is set by $\beta_{\mathrm{RL}}$, with $\beta_{\mathrm{RL}}=0$ corresponding to random actions and $\beta_{\mathrm{RL}}=\infty$ corresponding to always taking greedy actions with respect to the current estimate of the $Q$-function. Here we use an external ‘learning’ temperature scale, the inverse of which, $\beta_{\mathrm{RL}}$, is linearly ramped down as the number of episodes progresses. In the replay training stage, which is also $40$ episodes long, we replay the best-encountered protocol up to the given episode. Through this procedure, when the next exploratory training stage begins again, the agent is biased to do variations on top of the best-encountered protocol, effectively improving it, until it reaches a reasonably good fidelity.

Two learning curves of the RL agent are shown in Fig. 2. Notice the occurrence of suboptimal protocols even during later episodes due to the stochasticity of the exploration schedule. During every episode, the agent takes the best action (w.r.t. its current knowledge/experience) with a finite probability, or else a random action is chosen. This prevents the algorithm from immediately getting stuck in a high-infidelity (i.e. a bad) minimum. To guarantee convergence of the RL algorithm, the exploration probability is reduced as the number of episodes progresses (cf. discussion above). This becomes manifest in Fig. 2, where after many episodes the deviations from the good protocols decrease. In the end, the agent learns the best-encountered protocol as a result of using the replay schedule which speeds up learning (as can be seen by the bad shots becoming rarer with increasing the number of episodes). We show only these learned protocols in Fig. 1b and Fig. 3 of the Supplemental Material sup .

III Phases of Quantum Control

III.1 Single Qubit Manipulation

To benchmark the application of RL to physics problems, consider first a two-level system described by

[TABLE]

where $S^{\alpha}$, are the spin-$1/2$ operators. This Hamiltonian comprises both integrable many-body and non-interacting translational invariant systems, such as the transverse-field Ising model, graphene and topological insulators. The initial $|\psi_{i}\rangle$ and target $|\psi_{\ast}\rangle$ states are chosen as the ground states of (3) at $h_{x}=-2$ and $h_{x}=2$, respectively. We verified that the applicability of RL does not depend on this specific choice. Although there exists an analytical solution to solve for the optimal protocol in this case Hegerfeldt (2013), it does not generalize to non-integrable many-body systems. Thus, studying this problem using RL serves a two-fold purpose: (i) we benchmark the protocols obtained by the RL agent demonstrating that, even though RL is a completely model-free algorithm, it still finds the physically meaningful solutions by constructing a minimalistic effective model on-the-fly. The learning process is shown in Video 7; (ii) We reveal an important novel perspective on the complexity of quantum state manipulation which, as we show below, generalizes to many-particle systems. While experimental set-ups studying single-qubit physics can readily apply multiple control fields (e.g also control fields in the $y$-direction) in order to test RL on a non-trivial problem with a known solution, we restrict the discussion to a single control parameter.

For fixed total protocol duration $T$, the infidelity $h_{x}(t)\mapsto I_{h}(T)=1-F_{h}(T)$ represents a “potential landscape”, the global minimum of which corresponds to the optimal driving protocol. For bang-bang protocols, the problem of finding the optimal protocol becomes equivalent to finding the ground state configuration of a classical Ising model with complicated interactions Day et al. (2018). We map out the landscape of local infidelity minima $\{h_{x}^{\alpha}(t)\}_{\alpha=1}^{N_{\mathrm{real}}}$ using Stochastic Descent (SD), starting from random bang-bang protocol configurations sup . To study the correlations between the infidelity minima as a function of the total protocol duration $T$, we define the correlator $q(T)$, closely related to the Edwards-Anderson order parameter for the existence of spin glass order Castellani and Cavagna (2005); Hedges et al. (2009), as

[TABLE]

where $\overline{h_{x}}(t)=N_{\mathrm{real}}^{-1}\sum_{\alpha=1}^{N_{\mathrm{real}}}h_{x}^{\alpha}(t)$ is the sample-averaged protocol. If the minima $\{h_{x}^{\alpha}(t)\}_{\alpha=1}^{N_{\mathrm{real}}}$ are all uncorrelated, then $\overline{h_{x}}(t)\equiv 0$, and thus $q(T)=1$. On the other hand, if the infidelity landscape contains only one minimum, then $\overline{h_{x}}(t)\equiv h_{x}(t)$ and $q(T)=0$. The behaviour of $q(T)$, and the maximum fidelity $F_{h}(T)$ found using SD, together with a qualitative description of the corresponding infidelity landscapes are shown in Fig. 1.

The control problem for the constrained qubit exhibits three distinct control phases as a function of the protocol duration $T$. If $T$ is greater than the quantum speed limit $T_{\mathrm{QSL}}\approx 2.4$, one can construct infinitely many protocols which prepare the target state with unit fidelity, and the problem is in the controllable phase $III$, c.f. Fig. 1. The red line in Fig. 1b (iii) shows an optimal protocol of unit fidelity found by the agent, whose Bloch sphere representation can be seen in Video 3. In this phase, there is a proliferation of exactly degenerate, uncorrelated global infidelity minima, corresponding to protocols of unit fidelity, and the optimization task is easy.

At $T=T_{\mathrm{QSL}}$, the order parameter $q(T)$ exhibits a non-analyticity, and the system undergoes a continuous phase transition to a correlated phase $II$. For times smaller than $T_{\mathrm{QSL}}$ but greater than $T_{c}$, the degenerate minima of the infidelity landscape recede to form a correlated landscape with many non-degenerate local minima, as reflected by the finite value of the order parameter $0<q(T)<1$. As a consequence of this correlated phase, there no longer exists a protocol to prepare the target state with unit fidelity, since it is physically impossible to reach the target state while obeying all constraints. The infidelity minimization problem is non-convex, and determining the best achievable (i.e. optimal) fidelity [a.k.a. the global minimum] becomes difficult. Figure 1b (ii)shows the best bang-bang protocol found by our computer agent (see Video 2 and sup for protocols with quasi-continuous actions). This protocol has a remarkable feature: without any prior knowledge about the intermediate quantum state nor its Bloch sphere representation, the model-free RL agent discovers that it is advantageous to first bring the state to the equator – which is a geodesic – and then effectively turns off the control field $h_{x}(t)$, to enable the fastest possible precession about the $z$-axis 333Notice that the agent does not have control over the $z$ field.. After staying on the equator for as long as optimal, the agent rotates as fast as it can to bring the state as close as possible to the target, thus optimizing the final fidelity for the available protocol duration.

Decreasing the total protocol duration $T$ further, we find a second critical time $T_{c}\approx 0.6$. For $T<T_{c}$, $q(T)\equiv 0$ and the problem has a unique solution, suggesting that the infidelity landscape is convex. This overconstrained phase is labelled $I$ in the phase diagram (Fig. 1a). For $T<T_{c}$, there exists a unique optimal protocol, even though the achievable fidelity can be quite limited, see Fig. 1b (i) and Video 1. Since the state precession speed towards the equator depends on the maximum possible allowed field strength $h_{x}$, it follows that $T_{c}\to 0$ for $|h_{x}|\to\infty$.

*Relation to Counter-Diabatic and Fast-Forward Driving.—*Promising analytical approaches to state manipulation have recently been proposed, known as Shortcuts to Adiabaticity Demirplak and Rice (2003, 2005, 2008); Berry (2009); del Campo (2013); Deffner et al. (2014); Kolodrubetz et al. (2017); Sels and Polkovnikov (2017); Theisen et al. (2017); Petiziol et al. (2018). They include ideas such as (i) fast-forward (FF) driving, which comprises protocol that excite the system during the evolution at the expense of gaining speed, before taking away all excitations and reaching the target state with unit probability, and (ii) counter-diabatic (CD) driving, which ensures transitionless dynamics by turning on additional control fields. In general any FF protocol is related to a corresponding CD protocol. While for complex many-body systems, it is not possible to construct the mapping between FF and CD in general, the simplicity of the single-qubit setup (3) allows to use CD driving to find a FF protocol Bukov et al. (2018a). For an unbounded control field $h_{x}(t)$, the FF protocol at the quantum speed limit has three parts which can be understood intuitively on the Bloch sphere: (i) an instantaneous delta-function kick to bring the state to the equator, (ii) an intermediate stage where the control field is off, $h_{x}(t)\equiv 0$, which allows the state to precess along the equator, and (iii) a complementary delta kick to bring the state from the equator straight to the target Bukov et al. (2018a). Whenever the control field is bounded, $|h_{x}|\leq 4$, these delta kicks are broadened and take extra time, thus increasing $T_{\mathrm{QSL}}$. If the RL algorithm finds a unit-fidelity protocol, it is by definition a FF one. Comparing FF driving to the protocol found by our RL agent [cf. Fig. 1b, see also paragraphs above], we find indeed a remarkable similarity between the RL and FF protocols.

III.2 Many Coupled Qubits

The above results raise the natural question of how much more difficult state manipulation is in more complex quantum models. To this end, consider a closed chain of $L$ coupled qubits, which can be experimentally realized with superconducting qubits Barends et al. (2016), cold atoms Simon et al. (2011) and trapped ions Jurcevic et al. (2014):

[TABLE]

We set $g\!=\!1$ to avoid the anti-ferromagnet to paramagnet phase transition, and choose the paramagnetic ground states of Eq. (5) at fields $h_{x}=-2$ and $h_{x}=2$ for the initial and target state, respectively. We verified that the conclusions we draw below do not depend on the choice of initial and target states, provided they both belong to the paramagnetic phase. The details of the control field $h_{x}(t)$ are the same as in the single qubit case, and we use the many-body fidelity both the reward and the measure of performance. In this paper, we focus on $L>2$. The two-qubit optimization problem was shown to exhibit an additional symmetry-broken correlated phase, see Ref. Bukov et al. (2018b).

Figure 3 shows the phase diagram of the coupled qubits model. First, notice that while the overconstrained-to-glassy critical point $T_{c}$ survives, the quantum speed limit critical point $T_{\mathrm{QSL}}$ is (if existent at all) outside the short protocol-time range of interest. Thus, the glassy phase extends over to long and probably infinite protocol durations, which offers an alternative explanation for the difficulty of preparing many-body states with high fidelity. The glassy properties of this phase are analyzed extensively in Ref. Day et al. (2018). Second, observe that, even though unit fidelity is no longer achievable, there exist nearly optimal protocols with extremely high many-body fidelity 444We expect the many-body overlap between two states to be exponentially small in the system size $L$. at short protocol durations. This fact is striking because the Hilbert space of our system grows exponentially with $L$ and we are using only one control field to manipulate exponentially many degrees of freedom in a short time. Nonetheless, it has been demonstrated that two states very close to each other in, or of equal, fidelity can possess sufficiently different physical properties or be very far in terms of physical resources Benedetti et al. (2013); Bina et al. (2014); Mandarino et al. (2016); Bukov et al. (2018b). Hence, one should be cautious when using the fidelity as a measure for preparing many-body states and exploring other possible reward functions for training RL agents is an interesting avenue for future research.

Another remarkable characteristic of the optimal solution is that for the system sizes $L\geq 6$ both $q(T)$ and $-L^{-1}\log F_{h}(T)$ converge to their thermodynamic limit values with no visible finite-size corrections sup . This is likely related to the Lieb-Robinson bound for information propagation which suggests that information should spread over approximately $JT=4$ sites for the longest protocol durations considered.

IV Variational Theory for Nearly-Optimal Protocols

An additional feature of the optimal bang-bang solution found by the agent is that the entanglement entropy of the half system generated during the evolution always remains small, satisfying an area law sup . This implies that the system likely follows the ground state of some local, yet a-priori unknown effective Hamiltonian Vidmar et al. (2015). This emergent behavior motivated us to use the best protocols found by ML to construct simple variational protocols consisting of just a few bangs. Let us now demonstrate how this works by giving specific examples which, to our surprise, capture the essence of the phase diagram of quantum control both qualitatively and quantitatively.

IV.1 Single Qubit

By carefully studying the optimal driving protocols the RL agent finds in the case of the single qubit, we find a few important features. Focussing for the moment on bang-bang protocols, in the overconstrained and correlated phases [cf. Fig. 1b and Videos 1–3], we recognize an interesting pattern: for $T<T_{c}$, as we explained in Sec. III, there is only one minimum in the infidelity landscape, which dictates a particularly simple form for the bang-bang protocol – a single jump at half the total protocol duration $T/2$. On the other hand, for $T_{c}\leq T\leq T_{\mathrm{QSL}}$, there appears a sequence of multiple bangs around $T/2$, which grows with increasing the protocol duration $T$. By looking at the Bloch sphere representation, see Videos 1–3, we identify this as an attempt to turn off the $h_{x}$-field, once the state has been rotated to the equator. This trick allows for the instantaneous state to be moved in the direction of the target state in the shortest possible distance [i.e. along a geodesic].

Hence, it is suggestive to try out a three-pulse protocol as an ansatz for the optimal solution, see Fig. 4a: the first (positive) pulse of duration $\tau^{(1)}/2$ brings the state to the equator. Then the $h_{x}$-field is turned off for a time $\tilde{\tau}^{(1)}=T-\tau^{(1)}$, after which a negative pulse directs the state off the equator towards the target state. Since the initial value problem is time-reversal symmetric for our choice if initial and target states, the duration of the third pulse must be the same as that of the first one. We thus arrive at a variational protocol, parametrised by $\tau^{(1)}$, see Fig. 4a.

The optimal fidelity is thus approximated by the variational fidelity $\mathcal{F}_{h}(\tau^{(1)},T-\tau^{(1)})$ for the trial protocol [Fig. 4a], and can be evaluated analytically in a straightforward manner:

[TABLE]

However, since the exact expression is rather cumbersome, we choose not to show it explicitly. Optimizing the variational fidelity at a fixed protocol duration $T$, we solve the corresponding transcendental equation to find the extremal value $\tau^{(1)}_{\mathrm{best}}$, and the corresponding optimal variational fidelity $\mathcal{F}_{h}(T)$, shown in Fig. 4b-c. For times $T\leq T_{c}$, we find $\tau^{(1)}=T$ which corresponds to $\tilde{\tau}^{(1)}=0$, i.e. a single bang in the optimal protocol. The overconstrained-to-correlated phase transition at $T_{c}$ is marked by a non-analyticity at $\tau^{(1)}_{\mathrm{best}}(T_{c})=T_{c}\approx 0.618$. This is precisely the minimal time the agent can take, to bring the state to the equator of the Bloch sphere, and it depends on the value of the maximum magnetic field allowed [here $h_{\mathrm{max}}=4$]. Figure 4d shows that, in the overconstrained phase, the fidelity is optimised at the boundary of the variational domain, although $\mathcal{F}_{h}(\tau^{(1)},T-\tau^{(1)})$ is a highly nonlinear function of $\tau^{(1)}$ and $T$.

For $T_{c}\leq T\leq T_{\mathrm{QSL}}$, the time $\tau^{(1)}$ is kept fixed [the equator being the only geodesic for a rotation along the $\hat{z}$-axis of the Bloch sphere], while the second pulse time $\tilde{\tau}^{(1)}$ grows linearly, until the minimum time $T_{\mathrm{QSL}}\approx 2.415$ is eventually reached. The minimum time is characterised by a bifurcation in our effective variational theory, as the corresponding variational infidelity landscape develops two minima, see Fig. 4b,d. Past that protocol duration, our simplified ansatz is no longer valid, and the system is in the controllable phase. Furthermore, a sophisticated analytical argument based on optimal control theory can give exact expressions for $T_{c}$ and $T_{\mathrm{QSL}}$ Hegerfeldt (2013), in precise agreement with the values we obtained. The Bloch sphere representation of the variational protocols in Fig. 1b (dashed blue lines) for the single qubit are shown in Videos 4-6.

To summarize, for the single qubit example, the variational fidelity $\mathcal{F}_{h}(T)$ agrees nearly perfectly with the optimal fidelity $F_{h}(T)$ obtained using SD and Optimal Control, cf. Fig. 1a. We further demonstrate that our variational theory fully captures the physics of the two critical points $T_{c}$ and $T_{\mathrm{QSL}}$ sup . Interestingly, the variational solution for the single qubit problem coincides with the global minimum of the infidelity landscape all the way up to the quantum speed limit sup .

IV.2 Many Coupled Qubits

Let us also discuss the variational theory for the many-body system. Consider first the same one-parameter variational ansatz from Sec. IV.1, see Fig. 5a. Since the variational family is one-dimensional, we shall refer to this ansatz as the 1D variational theory. The dashed black line in Fig. 5c shows the corresponding 1D variational fidelity. We see that, once again, this ansatz captures correctly the critical point $T_{c}$ separating the overconstrained and the glassy phases. Nevertheless, a comparison with the optimal fidelity [see Fig. 5c] reveals that this variational ansatz breaks down in the glassy phase, although it rapidly converges to the optimal fidelity with decreasing $T$. Looking at Fig. 5b, we note that the value $\tau^{(1)}_{\mathrm{best}}$, which maximizes the variational fidelity, exhibits a few kinks. However, only the kink at $T=T_{c}$ captures a physical transition of the original control problem, while the others appear as artefacts of the simplified variational theory, as can be seen by the regions of agreement between the optimal and variational fidelities.

Inspired by the structure of the protocols found by our RL agent once again, see Video 8, we now extend the qubit variational protocol, as shown in Fig. 5d. In particular, we add two more pulses to the protocol, retaining its symmetry structure: $h_{x}(t)=-h_{x}(T-t)$, whose length is parametrised by a second, independent variational parameter $\tau^{(2)}/2$. Thus, the pulse length where the field is set to vanish, is now given by $\tilde{\tau}=T-\tau^{(1)}-\tau^{(2)}$. These pulses are reminiscent of spin-echo protocols, and appear to be important for entangling and disentangling the state during the evolution. Notice that this extended variational ansatz includes by definition the simpler ansatz from the single qubit problem discussed above, by setting $\tau^{(2)}=0$.

Let us now turn on the second variational parameter $\tau^{(2)}$, and consider the full two-dimensional variational problem:

[TABLE]

For the maximum-fidelity variational protocol, we show the best variational fidelity $\mathcal{F}^{\mathrm{2D}}_{h}$ [Fig. 5f] and the corresponding values of $\tau^{(1)}_{\mathrm{best}}$ and $\tau^{(2)}_{\mathrm{best}}$ [Fig. 5e]. There are two important points here: (i) Fig. 5f shows that the 2D variational fidelity seemingly reproduces the optimal fidelity on a much longer scale compared to the 1D variational ansatz, i.e. for all protocol durations $T\lesssim 3.3$. (ii) the 2D variational ansatz reduces to the 1D one in the overconstrained phase $T\leq T_{c}$. In particular, both pulse lengths $\tau^{(1)}_{\mathrm{best}}$ and $\tau^{(2)}_{\mathrm{best}}$ exhibit a non-analyticity at $T=T_{c}$, but also at $T^{\prime}\approx 2.5$. Interestingly, the 2D variational ansatz captures the optimal fidelity on both sides of $T^{\prime}$ which suggests that there is likely yet another transition within the glassy phase, hence the different shading in the many-body phase diagram [Fig. 3]. Similar to the 1D variational problem, here we also find artefact transitions [non-analytic behavior in $\tau^{(i)}_{\mathrm{max}}$ outside of the validity of the variational approximation].

In summary, in the many-body case, the same one-parameter variational ansatz only describes the behaviour in the overconstrained phase, cf. Fig. 3 (dashed line), up to and including the critical point $T_{c}$, but fails for $T>T_{c}$. Nevertheless, a slightly modified, two-parameter variational ansatz, motivated again by the solutions found by the ML agent (see Video 8), appears to be fully sufficient to capture the essential features of the optimal protocol much deeper into the glassy phase, as shown by the $\mathcal{F}^{\mathrm{2D}}_{h}(T)$ curve in Fig. 3. This many-body variational theory features an additional pulse, reminiscent of spin-echo, which appears to control and suppress the generation of entanglement entropy during the drive sup . Indeed, while the two-parameter ansatz is strictly better than the single-parameter protocol for all $T>T_{c}$, the difference between the two grows slowly as a function of time. It is only at a later time, $T\approx 1.3$, that the effect of the second pulse really kicks in, and we observe the largest entanglement in the system for the optimal protocol.

Using RL, we identified nearly-optimal control protocols Gingrich et al. (2016) which can be parametrized by a few d.o.f. Such simple protocols have been proven to exist in weakly-entangled one-dimensional spin chains Lloyd and Montangero (2014). However, the proof of the existence does not imply that these d.o.f. are easy to identify. Initially, the RL agent is completely ignorant about the problem and explores many different protocols, while it tries to learn the relevant features. In contrast, optimal control methods, such as CRAB Caneva et al. (2011), usually have a much more rigid framework, where the d.o.f. of the method are fixed from the beginning. This can limit the performance of those methods below the quantum speed limit Sørensen et al. (2016); sup .

One might wonder how the nearly-optimal protocols found using RL and SD correlate with the best variational protocols. For the problem under consideration, averaging parts of the set of bang-bang protocols, which contains randomly generated local minima of the infidelity landscape, $\{h_{x}^{\alpha}(t)\}_{\alpha=1}^{N_{\mathrm{real}}}$, (see insets in Fig. 7) results in protocols which resemble the continuous ones we found using GRAPE. The variational solutions are indeed close to these averaged solutions, although they are not exactly the same, since the variational protocols are constrained to take on three discrete values (positive, zero and negative), while the averaged protocols can take on any values in the interval $[-4,4]$. The RL agent cannot find these variational solutions because we have limited the actions space to having $h_{x}$ take the minimum or maximum allowable value and there is no way to take an action where $h_{x}=0$.

We also showed how, by carefully studying the driving protocols found by the RL agent, one can obtain ideas for effective theories which capture the essence of the underlying physics. This approach is similar to using an effective $\phi^{4}$-theory to describe the physics of the Ising phase transition. The key difference is that the present problem is out-of-equilibrium, where no general theory of statistical mechanics exists so far. We hope that, an underlying pattern between such effective theories can be revealed with time, which might help shape the guiding principles of a theory of statistical physics away from equilibrium.

V Glassy Behaviour

It is quite surprising that the dynamics of a non-integrable many-body quantum system, associated with the optimal protocol, is so efficiently captured by such a simple, two-parameter variational protocol, even in the regimes where there is no obvious small parameter and where spin-spin interactions play a significant role. Upon closer comparison of the variational and the optimal fidelities, one can find regions in the glassy phase where the simple variational protocol outperforms the numerical ‘best’ fidelity, cf. Fig. 3.

To better understand this behavior, we choose a grid of $N_{T}=28$ equally-spaced time steps, and compute all $2^{28}$ bang-bang protocols and their fidelities. The corresponding density of states (DOS) in fidelity space is shown in Fig. 6 for two choices of $T$ in the overconstrained and glassy phase. This allows us to unambiguously determine the ground state of the infidelity landscape (i.e. the optimal protocol). Starting from this ground state, we then construct all excitations generated by local in time flips of the bangs of the optimal protocol. The fidelity of the “1-flip” excitations is shown using red circles in Fig. 6. Notice how, in the glassy phase, these $28$ excitations have relatively low fidelities compared to the ground state, and are surrounded by $\sim 10^{6}$ other states. This has profound consequences: as we are ‘cooling’ down in the glassy phase, searching for the optimal protocol and coming from a state high up in the infidelity landscape, if we miss one of the $28$ elementary excitations, it becomes virtually impossible to reach the global ground state and the situation becomes much worse if we increase the number of steps $N_{T}$. On the contrary, in the overconstrained phase, the smaller value of the DOS at the “1-flip” excitation ($\sim 10^{2}$) makes it easier to reach the ground state.

The green crosses in Fig. 6 show the fidelity of the “2-flip” excitations. By the above argument, a “2-flip” algorithm would not see the phase as a glass for $T\lesssim 2.5$, yet it does so for $T\gtrsim 2.5$, marked by the different shading in Fig. 3. Correlated with this observation, we find a signature of a transition also in the improved two-parameter variational theory in the glassy phase [see Sec. IV.2 and kinks at $T^{\prime}$ in Fig. 5e]. In general, we expect the glassy phase to exhibit a series of phase transitions, reminiscent of the random $k$-SAT problems Mézard et al. (2002); Morampudi et al. (2017). The glassy nature of this correlated phase has been studied in detail in Ref.Day et al. (2018), by mapping this optimal control problem to an effective classical spin energy function which governs the control phase transitions.

In contrast to the single-qubit system, there are also multiple attractors present in the glassy phase of the many-body system [see Fig. 7]. Each attractor has a typical representative protocol [Fig. 7 insets]. Even though intra-attractor protocols share the same averaged profile, they can nevertheless have a small mutual overlap – comparable to the overlap of inter-attractor protocols. This indicates that in order to move in between protocols within an attractor, highly non-local moves are necessary. For this reason, GRAPE Khaneja et al. (2005b), an algorithm which performs global updates on the protocol by computing exact gradients in the control landscape, also performs very well on our optimisation problem. Similar to SD, in the glassy phase GRAPE cannot escape local minima in the infidelity landscape and, therefore, the same three attractors are found with comparable relative populations to SD, but intra-attractor fluctuations are significantly suppressed due to GRAPE’s non-local character.

VI Outlook & Discussion

In this work we demonstrated the usefulness of Q-Learning to manipulate single-particle and many-body quantum systems. Q-Learning is only one of many Reinforcement Learning algorithms, including SARSA, Policy Gradient and Actor Critic methods, just to name a few. In the Supplemental Material we showed that Q-Learning’s performance is comparable to many of the leading Optimal Control algorithms sup . It will be interesting and desirable to compare different RL algorithms among themselves on physical quantum systems. An exciting future direction is to investigate which advantages Deep Learning offers in the context of quantum control, and there exist recent studies exploring Deep RL in a physics Fösel et al. (2018); August and Hernández-Lobato (2018); Niu et al. (2018).

Looking forward to controlling non-integrable many-body systems, an important question arises as to how the computational cost of Q-Learning scales with the system size $L$. As we explained in Sec. II, the Q-Learning algorithm can be decomposed into a ‘learning’ part, and an ‘interaction with the environment’ part where all physics/dynamics happens. The learning part does not know about the state of the quantum system – it only keeps track of the value of the magnetic field at a given time $(t,h_{x}(t))$. As a result of this choice, for a single global drive the ‘learning part’ of the algorithm is independent of the system size $L$ since it depends only on a single scalar reward - the fidelity of the final state. The RL algorithm is instead computationally limited by the size of the action and state spaces. As currently implemented, this means that the RL algorithm is limited to finding short protocols (since the state space scales exponentially with the number of bangs). However, it may be possible to circumvent this bottleneck by using Deep RL which uses neural networks to represent the Q-function.

One place where the system size implicitly enters the computational costs of the RL protocol is through the number of episodes needed to train the RL algorithm. At every time step, one solves Schrödinger’s equation to simulate the dynamics. The solver’s scaling with $L$ depends on how the time evolution is implemented: in spin-$1/2$ systems, for exact diagonalization [used here] the computational cost scales exponentially $2^{2L}$, while a more sophisticated Krylov method alleviates this somewhat to $L^{2}2^{L}$ Saad (1992), and matrix product states only scale as $L^{2}$ [in the worst case] García-Ripoll (2006). Therefore, combining RL with existing approximate techniques to evolve quantum states can lead to a significant reduction of CPU time, provided applying these techniques is justified by the underlying physics.

The present work demonstrates the suitability of RL for manipulating/controlling quantum systems. Yet, it does not explore how one can improve the Q-Learning algorithm, and adjust it to the specific needs of quantum control. Let us briefly list a few possible directions that the interested reader may want to keep in mind: (i) Alternative definitions of the RL state space, see Ref. Sutton and Barto (2017), may prove advantageous, depending on the need of the problem setup, since this defines the agent’s knowledge about the physical system. For instance, if the RL agent is to be coupled to an experiment, one cannot use the wavefunction for this purpose, whereas wavefunctions may be accessible in numerical simulations. We find that the choice of RL state space influences the learning capabilities of the RL agent.(ii) Another way to increase performance is to add more controls. This only increases the possibility to reach a higher fidelity, but it comes at a cost of a potential slow-down, due to a higher computational demand to explore the increased RL state space. (iii) In addition, choosing a suitable family of protocols and how to parametrize it, may also lead to increased performance in RL. We used bang-bang protocols because of their computational simplicity, yet the needs of a given problem may justify another choice: the experimental realization of bang-bang protocols is limited by the resolution with which a pulse can be stabilized, which is set by the experimental apparatus. Alternatively, the RL setup can be formulated to control the size of some generalized Fourier coefficients, an idea underlying the CRAB algorithm. (iv) On the algorithmic side, one can also optimize the exploration and replay schedules which control the learning efficiency with increasing the number of training episodes, and influence the RL agent’s learning speed.

Reinforcement Learning algorithms are versatile enough and can be suitably combined with existing ones. For instance, applying RL to complex problems with glassy landscapes is likely to benefit from a pre-training stage. Such a beneficial behaviour has already been observed in the context of deep RL Mnih et al. (2015); Silver et al. (2016). For the purpose of pre-training, in certain cases it may be advantageous to combine RL with existing derivative-based optimal control methods, such as GRAPE and CRAB, or even exhaustive search, so that one starts the optimization from a reasonable ‘educated guess’. In the recent years, it was shown that derivative-based and feedback-loop control methods can be efficiently combined to boost performance Egger and Wilhelm (2014). Vice-versa, RL’s exploration schedule defined on a suitable abstract RL-state space, may prove a useful addition to improve on already existing algorithms.

Using RL, we revealed the existence of control phase transitions, and showed their universality in the sense that they also affect the behaviour of state-of-the-art optimal control methods. The appearance of a glassy phase, which dominates the many-body physics, in the space of protocols of the quantum state manipulation problem, could have far-reaching consequences for efficiently manipulating systems in condensed matter experiments. Quantum computing relies heavily on our ability to prepare states with high fidelity, yet finding high efficiency state manipulation routines remains a difficult problem. Highly controllable quantum emulators, such as ultracold atoms and ions, depend almost entirely on the feasibility to reach the correct target state, before it can be studied. We demonstrated how, a model-free RL agent can provide valuable insights in constructing variational theories which capture almost all relevant features of the dynamics generated by the optimal protocol. Unlike the optimal bang-bang protocol, the simpler variational protocol is robust to small perturbations, while giving comparable fidelities. This implies the existence of nearly optimal protocols, which do not suffer from the exponential complexity of finding the global minimum of the entire optimization landscape. Finally, in contrast with optimal control methods such as SGD, GRAPE, and CRAB that assume an exact model of the physical system, the model-free nature of RL suggests that it can be used to design protocols even when our knowledge of the physical systems we wish to control is incomplete or our system is noisy or disordered Bukov (2018).

The existence of phase transitions in quantum control problems may have profound consequences beyond physical systems. We suspect that the glassy behavior observed here maybe a generic feature of many control problems and it will be interesting to see if this is indeed the case. It is our hope that given the close connections between optimal control and RL, the physical interpretation of optimization problems in terms of a glassy phase will help in developing novel efficient algorithms and help spur new ideas in RL and artificial intelligence.

Acknowledgements.

We thank J. Garrahan, M. Heyl, M. Schiró and D. Schuster for illuminating discussions. MB, PW and AP were supported by NSF DMR-1813499, ARO W911NF1410540 and AFOSR FA9550-16-1-0334. AD is supported by a NSERC PGS D. AD and PM acknowledge support from Simon’s Foundation through the MMLS Fellow program. DS acknowledges support from the FWO as post-doctoral fellow of the Research Foundation – Flanders and CMTV. We used QuSpin for simulating the dynamics of the qubit systems Weinberg and Bukov (2017, 2018). The authors are pleased to acknowledge that the computational work reported on in this paper was performed on the Shared Computing Cluster which is administered by Boston University’s Research Computing Services. The authors also acknowledge the Research Computing Services group for providing consulting support which has contributed to the results reported within this paper.

S1 Glassy Behaviour of Different Machine Learning and Optimal Control Algorithms with Local and Nonlocal Flip Updates

S1.1 Stochastic Descent

To benchmark the results obtained using Reinforcement Learning (RL), we use a greedy stochastic descent (SD) algorithm to sample the infidelity landscape minima containing the driving protocols. We restrict our SD algorithm to exploring bang-bang protocols, for which $h_{x}(t)\in\{\pm 4\}$. The algorithm starts from a random protocol configuration and proposes local field updates at a time $t$ chosen uniformly in the interval $[0,T]$. The updates consist in changing the applied field $h_{x}(t)\rightarrow h_{x}^{\prime}(t)$ only if this increases the fidelity. Ideally, the protocol is updated until all possible local field updates can only decrease the fidelity. Practically, for some protocol durations, ensuring that a true local minima with respect to $1$-flip is reached can be computational expensive. Therefore, we restrict the number of fidelity evaluations to be at most $20\times T/\delta t$. In this regard, the obtained protocol is a local minimum with respect to local ($1$-flip) field updates. The stochastic descent is repeated multiple times with different initial random protocols. The set of protocols $\{h^{\alpha}|\alpha=1,\dots,N_{\mathrm{real}}\}$ obtained with stochastic descent is used to calculate the glass-like order parameter $q(T)$ (see main text). A Python implementation of the algorithm is available on Github.

S1.2 CRAB

Chopped RAndom Basis (CRAB) is a state-of-the-art optimal control algorithm designed to tackle many-body quantum systems Doria et al. (2011); Caneva et al. (2011). The idea behind CRAB is to decompose the unknown driving protocol into a complete basis (Fourier, Laguerre, etc.), and impose a cut-off on the number of ‘harmonics’ kept for the optimisation. The algorithm then uses an optimiser to find the values for the expansion coefficients, which optimise the cost function of the problem.

Following Ref. Caneva et al. (2011), we make a Fourier-basis ansatz for the driving protocol.

[TABLE]

where the Fourier coefficients $\{A_{i},B_{i},\omega_{i}\}$ which parametrise the protocol are found using an implementation of the Nelder-Mead optimization method in the SciPy python library. The number of harmonics kept in the optimisation is given by $N_{c}$. The CRAB algorithm uses two auxiliary functions, defined by the user: the first function $h_{0}(t)$ is a trial initial guess ansatz for the protocol, while the second function, $\lambda(t)$ imposes the boundary conditions $\lambda\to\infty$ for $t\to 0$ and $t\to T$ to the Fourier expansion term.

The cost function which we optimise in the state manipulation problem

[TABLE]

contains the fidelity $\mathcal{F}(\{A_{i},B_{i},\omega_{i}\})$ at the end of the protocol, and an additional penalty coming from the $L^{2}$ norm of the protocol to keep the optimal protocols bounded. The last constraint is required for a better and honest comparison with the RL, SD, and GRAPE algorithms.

Applying CRAB to the state manipulation problem from the main text systematically, we choose $\lambda(t)=1/\sin(\pi t/T)^{2}$, and $h_{0}(t)=-2+4t/T$. We also consider $N_{c}=10,20$ to study how much an effect increasing the number of degrees of freedom will have on the optimal protocols found. For each value of $N_{c}$ we start the optimization algorithm with 10 random initial configurations. We define the optimal protocol as the protocol with the best fidelity out of that group of 10. The random initial frequencies $\omega_{i}$ are chosen the same way as outlined in Ref. Caneva et al. (2011) while the amplitudes $A_{i}$ and $B_{i}$ are chosen uniformly between $-10$ and $10$.

S1.3 GRAPE

GRadient Ascend Pulse Engineering, is a numeric derivative-based optimal control method, first introduced in the context of NMR spectroscopy Khaneja et al. (2005b). As suggested by its name, the method performs gradient optimization. Instead of restricting the protocols to bang-bang type, the method works with quasi-continuous protocols. Protocol magnitudes can take on any value within the allowed manifold but, unlike CRAB, are piecewise constant in time.

In the present case, one can efficiently compute the gradient of the fidelity as follows. Consider the fidelity for some trial protocol $h_{x}(t)$,

[TABLE]

where $U(T,0)$ denotes the time evolution operator from [math] to $T$. Let us further decompose the Hamiltonian as $H=H_{0}+h_{x}(t)X$, where $H_{0}$ is the part over which we have no control, and $X$ denotes the operator we control. The functional derivative of the fidelity with respect to the protocol thus becomes

[TABLE]

Although this expression appears hard to evaluate, it takes on a very simple form

[TABLE]

where $\left|\psi_{i}(t)\right\rangle=U(t,0)\left|\psi_{i}\right\rangle$ denotes the initial state propagated to time $t$ and $\left\langle\phi(t)\right|=\left\langle\psi_{i}(T)|\psi_{\ast}\right\rangle\left\langle\psi_{\ast}\right|U(T,t)$ denotes the (scaled) final state propagated back in time to time $t$. Notice that this procedure requires us to exactly know this time evolution operator (i.e. a model of physical system to be controlled). Hence, by propagating both the initial state forward and the target state backward in time one gets access to the full gradient of the control landscape.

To find a local maximum one can simply now gradient ascend the fidelity. A basic algorithm thus goes as follows:

1. (i)

Pick a random initial magnetic field $h^{0}_{x}(t)$. 2. (ii)

Compute first $\left|\psi_{i}(t)\right\rangle$ and then $\left\langle\phi(t)\right|$ for the current setting of the magnetic field $h^{N}_{x}(t)$. 3. (iii)

Update the control field $h^{N+1}_{x}(t)=h^{N}_{x}(t)+\epsilon_{N}{\rm Im}\left[\left\langle\phi(t)\right|X\left|\psi_{i}(t)\right\rangle\right]$. Note that this step can be upgraded to a second-order Newton method to improve the performance of the algorithm, see Refs. Machnes et al. (2011); Agundez et al. (2017); De Fouquieres et al. (2011). 4. (iv)

Repeat (ii) and (iii) until the desired tolerance is reached.

Here $\epsilon_{N}$ is the step size in each iteration. Choosing a proper step size can be a difficult task. In principle the fidelity should go up after each iteration but if the step size it too large, the algorithm can overshoot the maximum, resulting in a worse fidelity. To avoid this, one should adapt the step size during the algorithm. Numerically we have observed that, in order to avoid overshooting saddles/maxima, $\epsilon_{N}\propto 1/\sqrt{N}$.

S1.4 Comparison between the RL, SD, CRAB and GRAPE

While a detailed comparison between the ML algorithms and other optimal control algorithms is an interesting topic, it is beyond the scope of the present paper. Below, we only show that both the RL algorithm is capable of finding optimal bang-bang protocols for the quantum control problem from the main text, and that its performance rivals that of SD and the state-of-the-art algorithms for many-body quantum problems CRAB and GRAPE.

The result for the single qubit $L=1$ are shown in Fig. S2 (left). The most important points can be summarised as follows:

• •

RL, GS and GRAPE all find the optimal protocol.

• •

Below the quantum speed limit, $T_{\mathrm{QSL}}$, CRAB finds good, but clearly suboptimal protocols. The plots also show the glassiness represents a generic feature of the constrained optimization problem and not the method used to perform the optimization. Increasing the cutoff $N_{c}$, and with it the number of effective degrees of freedom, does not lead to a sizeable improvement in CRAB. One explanation is the following: for the single qubit, we know that the variational protocol, which contains at most two bangs (see Fig. 1b), is a global minimum of the optimisation landscape. Such a protocol can easily be approximated using up to $N_{c}=20$ harmonics. However, allowing more degrees of freedom comes at a huge cost due to the glassiness of the problem: there exist many quasi-degenerate local minima for the algorithm to get stuck in.

The comparison for the Many Coupled Qubits system for $L=6$ are shown in Fig. S2 (right). Larger values of $L$ do not introduce any change in the behaviour, as we argue in a subsequent section below. In the many-body case, the variational protocol, which contains four bangs, is shown not to be the global minimum of the infidelity landscape. Instead, the true global minimum contains many more bangs which, however, only marginally improve the fidelity.

• •

All algorithms give reasonable fidelities, see Fig. S2 (right)

• •

Even though GRAPE seems to display better performance out of the four methods, one should not forget that this algorithm, which uses global flips, requires knowledge of all fidelity gradients – valuable information which is not easily accessible through experimental measurements. One has to keep in mind though, that this comparison is not completely honest, since GRAPE allows for the control field $h_{x}$ to take any value in the interval $[-4,4]$, which offers a further advantage over the bang-bang based RL. On the other hand, the model-free RL rivals the performance of GRAPE at all protocol durations, and outperforms CRABS, even for the single qubit below the quantum speed limit, where the problem enters the glassy phase. The reason for the seemingly slowly decreasing performance of RL with $T$ is that, since the total protocol duration is held fixed, the number of bangs increases exponentially with $T\sim N_{T}$, and hence the state space which has to be explored by the agent also grows exponentially. A detailed study of the scaling is postponed to future studies.

• •

All algorithms suffer from the glassiness in the optimisation landscape. This is not surprising, since the glass phase is an intrinsic property of the infidelity landscape, as defined by the optimisation problem, and does not depend on which algorithm is used to look for the optimal solution. Fig. S1 shows the fidelity traces as a function of running time for GRAPE. Comparing this to the corresponding results for SD, see Fig. 7, we see a strikingly similar behaviour, even though GRAPE uses nonlocal flip updates in contrast to SD. This means that, in the glassy phase, GRAPE also gets stuck in suboptimal attractors, similar to SD and RL. Thus, as an important consequence, the glassy phase affects both local and nonlocal-update algorithms.

S2 Performance of the Different Driving Protocols for the Qubit

It is interesting to compare the bang-bang and quasi-continuous driving protocols found by the agent to a simple linear protocol, which we refer to as Landau-Zener (LZ), and the geodesic protocol, which optimizes local fidelity close to the adiabatic limit essentially slowing down near the minimum gap Tomka et al. (2016). We find that the RL agent offers significantly better solutions in the overconstrained and glassy phases, where the optimal fidelity is always smaller than unity. The Hamiltonian of the qubit together with the initial and target states read:

[TABLE]

where $|\psi_{i}\rangle$ and $|\psi_{\ast}\rangle$ are the ground state of $H(t)$ for $h_{i}=-2$ and $h_{\ast}=+2$ respectively. Note that for bang-bang protocols, the initial and target states are not eigenstates of the control Hamiltonian since $h_{x}(t)$ takes on the values $\pm 4$.

The RL agent is initiated at the field $h(t=0)=h_{\mathrm{min}}=-4.0$. The RL protocols are constructed from the following set of jumps, $\delta h_{x}$, allowed at each protocol time step $\delta t$:

• •

bang-bang protocol: $\delta h_{x}\in\{0.0,\pm 8.0\}$ which, together with the initial condition, constrains the field to take the values $h_{x}(t)\in\{\pm 4.0\}$.

• •

quasi-continuous protocol: $\delta h_{x}\in\{0.0,\pm 0.1,\pm 0.2,\pm 0.5,\pm 1.0,\pm 2.0,\pm 4.0,\pm 8.0\}$. We restrict the actions available in a state to ensure $h_{x}(t)\in[-4.0,4.0]$.

Interestingly, the RL agent figures out that it is always advantageous to first jump to $h_{\mathrm{max}}=+4.0$ before starting the evolution, as a consequence of the positive value of the coefficient in front of $S^{z}$.

The analytical adiabatic protocols are required to start and end in the initial and target states, which coincide with the ground states of the Hamiltonians with fields $h_{i}=-2.0$ and $h_{\ast}=2.0$, respectively. They are defined as follows:

• •

Landau-Zener(LZ) protocol: $h_{x}(t)=(h_{\ast}-h_{i})t/T+h_{i}$

• •

geodesic protocol: $h_{x}(t)=\tan(at+b)$, where $b=\arctan(h_{i})$ and $a=\arctan(h_{\ast}-b)/T$.

Figure S3 shows a comparison between these four protocol types for different values of $T$, corresponding to the three quantum control phases. Due to the instantaneous gap remaining small compared to the total protocol duration, the LZ and geodesic protocols are very similar, irrespective of $T$. The two protocols significantly differ only at large $T$, where the geodesic protocol significantly outperforms the linear one. An interesting general feature for the short protocol durations considered is that the fidelities obtained by the LZ and geodesic protocols are clearly worse than the ones found by the RL agent. This points out the far-from-optimal character of these two approaches, which essentially reward staying close to the instantaneous ground state during time evolution. Looking at the fidelity curves in Fig. S3, we note that, before reaching the optimal fidelity at the end of the ramp for the overconstrained and glassy phases, the instantaneous fidelity drops below its initial value at intermediate times. This suggests that the angle between the initial and target states on the Bloch sphere becomes larger in the process of evolution, before it can be reduced again. Such situation is very reminiscent of counter-diabatic or fast forward driving protocols, where the system can significantly deviate from the instantaneous ground state at intermediate times Sels and Polkovnikov (2017); Zhou et al. (2017); Jarzynski et al. (2017). Such problems, where the RL agent learns to sacrifice local rewards in view of obtaining a better total reward in the end are of particular interest in RL Sutton and Barto (2017).

S3 Critical Scaling Analysis of the Control Phase Transitions

In this section, we show convincing evidence for the existence of the phase transitions in the quantum state manipulation problem discussed in the main text. We already argued that there is a phase transition in the optimization problem as a function of the protocol time $T$. Mathematically the problem is formulated as

[TABLE]

i.e. as finding the optimal driving protocol $h^{\mathrm{optimal}}_{x}(t)$, which transforms the initial state (the ground state of the Hamiltonian $H$ at $h_{x}=-2$) into the target state (the ground state of $H$ corresponding to $h_{x}=2$) maximizing the fidelity in protocol duration $T$ under unitary Schrödinger evolution. We assume that the ramp protocol is bounded, $h_{x}(t)\in[-4,4]$, for all times during the ramp. In this section, we restrict the analysis to bang-bang protocols only, for which $h_{x}(t)\in\{\pm 4\}$. The minimum protocol time step is denoted by $\delta t$. There are two different scaling limits in the problem. We define a continuum limit for the problem as $\delta t\to 0$ while keeping the total protocol duration $T=const$. Additionally, there is the conventional thermodynamic limit, where we send the system size $L\to\infty$.

As we already alluded to in the main text, one can think of this optimization problem as a minimization in the infidelity landscape, determined by the mapping $h_{x}(t)\mapsto I_{h}(T)=1-F_{h}(T)$, where each protocol is assigned a point in fidelity space – the probability of being in the target state after evolution for a fixed protocol duration $T$. Finding the global minimum of the landscape then corresponds to obtaining the optimal driving protocol for any fixed $T$.

To obtain the set of local minima $\{h_{x}^{\alpha}(t)|\alpha=1,\dots,N_{\mathrm{real}}\}$ of the infidelity landscape at a fixed total protocol duration $T$ and protocol step size $\delta t$, we apply Stochastic Descent(SD), see above, starting from a random protocol configuration, and introduce random local changes to the bang-bang protocol shape until the fidelity can no longer be improved. This method is guaranteed to find a set of representative local infidelity minima with respect to “$1$-flip” dynamics, mapping out the bottom of the landscape of $I_{h}(T)$. Keeping track of the mean number of fidelity evaluations $N_{\mathrm{fidevals}}$ required for this procedure, we obtain a measure for the average time it takes the SD algorithm to settle in a local minimum. While the order parameter $q(T)$ (see below) was used in the main text as a measure for the static properties of the infidelity landscape, dynamic features are revealed by studying the number of fidelity evaluations $N_{\mathrm{fidevals}}$.

As discussed in the main text, the rich phase diagram of the problem can also be studied by looking at the order parameter function $q$ (closely related to the Edwards-Anderson order parameter for detecting glassy order in spin systems Castellani and Cavagna (2005)):

[TABLE]

Here, $N_{T}$ is the total number of protocol time steps of fixed width $\delta t$, $N_{\mathrm{real}}$ is the total number of random protocol realisations $h_{x}^{\alpha}(t)$ probing the minima of the infidelity landscape (see previous paragraph), and the factor $1/16$ serves to normalise the squared bang-bang drive protocol $h_{x}^{2}(t)$ within the range $[-1,1]$.

S3.1 Single Qubit

For $T>T_{\mathrm{QSL}}$, the optimization problem of the single qubit ($L=1$) is controllable, and there exist infinitely many protocols which can prepare the target state with unit fidelity. In analogy with the random $k$-SAT problem Mézard et al. (2002), we call this the controllable (or underconstrained) phase of the quantum control problem. Intuitively, this comes about due to the large total protocol durations available which allow one to correct a set of ‘wrong moves’ at a later time in order to achieve a better fidelity at the end of the ramp. We have seen that both the Reinforcement Learning (RL) and Stochastic Descent (SD) agents readily and quickly learn to exploit this feature for optimal performance. In this phase, which is not present in the thermodynamic limit $L\to\infty$, there exist analytical results to compute driving protocols of unit fidelity based on the geodesic and counter-diabatic approaches Hegerfeldt (2013); Kolodrubetz et al. (2017); Sels and Polkovnikov (2017). The driving protocols $h_{x}^{\alpha}(t)$, corresponding to the minima of the infidelity landscape $I_{h}(T)$, are completely uncorrelated, resulting in $\overline{h_{x}}(t)=0$ and, thus, $q=1$. As $T\searrow T_{\mathrm{QSL}}$, the infidelity minima start becoming correlated, reflected in a drop in the value of $q(T)$. At the time $T=T_{\mathrm{QSL}}$, a phase transition occurs to a glassy phase with shallow, quasi-degenerate infidelity minima corresponding to many almost equally optimal protocols. Fig. S4a shows that the order parameter $q(T)$ develops a clear non-analyticity in the continuum limit $\delta t\to 0$, which proves the existence of a phase transition in protocol space. At the same critical time $T_{\mathrm{QSL}}$, a drastic rapid increase is detected in the number of fidelity evaluations required to map out the infidelity minima, see Fig. S4c.

For $T_{c}<T<T_{\mathrm{QSL}}$ the control problem is in the glassy phase. We showed in the main text by examining the DOS of protocols with respect to local flips of the field, that finding the optimal protocol appears as complicated as finding the ground state of a glass. This is reflected in the observed increase of the average number of fidelity evaluations $N_{\mathrm{fidevals}}$ with decreasing $\delta t$ (c.f Fig. S4c), and a decrease in the learning rate of the RL and SD agents. The local minima protocols $\{h_{x}^{\alpha}(t)\}$ are strongly correlated, as can be seen from the finite value of the order parameter $q(T)$ in Fig. S5a. More importantly, for any practical purposes, unit fidelity can no longer be obtained under the given dynamical constraints.

When we reach the second critical point $T=T_{c}$, another phase transition occurs from the glassy to an overconstrained phase. At $T=T_{c}$, the order parameter reaches zero, suggesting that the infidelity landscape contains a single minimum. In this phase, i.e. for $T<T_{c}$, the protocol duration is too short to achieve a good fidelity. Nonetheless, in the continuum limit $\delta t\to 0$, there exists a single optimal protocol, although the corresponding maximum fidelity is far from unity. In this overconstrained phase, the optimization problem becomes convex and easy to solve. This is reflected by the observation both the optimal quasi-continuous and bang-band protocols found by the RL agent are nearly identical, cf. Fig. S3. The dynamic character of the phase transition is revealed by a sudden drop in the number of fidelity evaluations $N_{\mathrm{fidevals}}$.

S3.2 Coupled Qubits

One can also ask the question what happens to the quantum control phases in the thermodynamic limit, $L\to\infty$. To this end, we repeat the calculation for a series of chain lengths $L$. We omit the case $L=2$, in which the physics of the control problem has a different character, exhibiting spontaneous symmetry breaking in the glasy phase Bukov et al. (2018b). Due to the non-integrable character of the many-body problem, we are limited to small system sizes. However, Fig. S5 shows convincing data that we capture the behaviour of the system in the limit $L\to\infty$ for the relatively short protocol durations under consideration. Moreover to our surprise the finite-size effects almost entirely disappear for $L\geq 6$ for all range of protocol durations we are considering. It seems that system is able to find an optimal solution, where the information simply does not propagate outside of a very small region and hence the optimal protocol rapidly becomes completely insensitive to the system size.

Figure S5c-d shows the system size scaling of the negative logarithmic many-body fidelity. While at $L=4$ we do see remnants of finite-size effects, starting from $L=6$ the curves are barely changing. A similar rapid system-size convergence is observed also for the order-parameter $q(T)$ (see Fig. S5a-b) and the entanglement entropy of the half chain (Fig. S5e). The protocol step size dependence of the order parameter $q(T)$, the average fidelity evaluations $N_{\mathrm{fidevals}}$, and the entanglement entropy $S_{\mathrm{ent}}^{L/2}$ of the half-chain are shown in Figs. S4b, S4d and S5f.

S4 Video Material

This paper is accompanied by eight short videos, labeled Video $1$ through $8$. The videos are available as part of the Supplemental Material published together with the paper. A legend for the videos can be found on https://mgbukov.github.io/RL_movies/.