QAOA implementation in a quantum reinforcement learning algorithm for CERN beam lines

QAOA implementation in a quantum reinforcement learning algorithm for CERN beam lines Master’s Thesis Master’s degree in Advanced Telecommunication Technologies (MATT) Author: Sebastià Nicolau Orell (3) Co-directors: Sofia Vallecorsa (1), Elías F. Combarro (2) &Eduard Alarcón (3) (1) CERN, 1, Esplanade des Particules, Geneva, CH 1211 (2) Computer Science Department, University of Oviedo, C. San Francisco, 3, 33003 Oviedo, Asturias, Spain (3) Escola Tècnica Superior d’Enginyeria de Telecomunicació de Barcelona, Universitat Politècnica de Catalunya, Carrer de Jordi Girona, 1-3, 08034 Barcelona, Spain Contents Acknowledgements ............................... i Abstract..................................... iii Listoffigures.................................. v 0 Introduction 1 1 Classical Reinforcement Learning 3 1.1 Fundamentals............................... 3 1.2 Q-learning................................. 4 1.3 DeepQ-learning.............................. 6 2 Quantum computing 7 2.1 Gate-based model fundamentals . . . . . . . . . . . . . . . . . . . . . 7 2.2 Adiabatic quantum computing . . . . . . . . . . . . . . . . . . . . . . 10 3 Quantum Machine Learning 13 3.1 Quantum Boltzmann Machine . . . . . . . . . . . . . . . . . . . . . . 13 3.2 Quantum Approximate Optimization Algorithm (QAOA) . . . . . . . 15 3.2.1 RecursiveQAOA ......................... 16 3.3 Quantum Reinforcement Learning . . . . . . . . . . . . . . . . . . . . 17 3.3.1 Energy-based Reinforcement Learning . . . . . . . . . . . . . 17 4 Environments 19 4.1 1D steering problem at CERN . . . . . . . . . . . . . . . . . . . . . . 19 3 5 Implementation 21 5.1 FERL algorithm implementation . . . . . . . . . . . . . . . . . . . . 21 5.2 QAOA implementation as solver . . . . . . . . . . . . . . . . . . . . . 24 6 Results 27 7 Conclusion 31 CONTENTS i Acknowledgements I would like to specially thank Elías Combarro for his constant support, Sofia Vallecorsa for giving me this opportunity, and Michael Schenk for guiding me through the implementation. I also want to thank Luis Meijueiro, Raúl Alonso and Andrés García from CTIC for their technical support. Finally, I would like to thank Eduard Alarcón for his accompaniement throughout this process. ii CONTENTS CONTENTS iii Abstract Quantum reinforcement learning (QRL), a cutting-edge field at the intersection of quantum computing and artificial intelligence, has the potential to revolutionize domains like chemistry and autonomous systems with its computational advantages. The free energy-based reinforcement learning algorithm (FERL) is a QRL algorithm inspired by the classical deep Q-learning algorithm (DQN) which uses a quantum Boltzmann machine (QBM) as the Q-function with quantum annealing (QA) as the solver. It has demonstrated considerable advantage compared to its classical counterpart DQN. The aim of this project is exploring the replacement of QA with the quantum approximate optimization algorithm (QAOA) and its recursive variant (RQAOA). We focus on the one-dimensional beam target steering control task based on the beam optics of the TT24-T4 transfer line at CERN. We perform and compare simulations of QA and QAOA, obtaining successful results for both QAOA and RQAOA, noticing that the performance of FERL using RQAOA for p=1 as a solver is very similar to the obtained through simulated QA for this scenario. The success of this experiment might lead the application of the QAOA solver to the hybrid actor-critic algorithm in development at CERN which allows to tackle control tasks using a continuous action space. iv CONTENTS List of Figures 1.1 Reinfocement learning scheme. . . . . . . . . . . . . . . . . . . . . . . 3 2.1 Bloch sphere; a geometrical representation of a two-level quantum system.................................... 8 2.2 Hadamardgate............................... 8 2.3 CNOTgate................................. 9 2.4 Quantum teleportation circuit. . . . . . . . . . . . . . . . . . . . . . . 10 3.1 Quantum machine learning diagram. . . . . . . . . . . . . . . . . . . 13 3.2 Clamped Boltzmann Machine used in FERL. . . . . . . . . . . . . . . 14 3.3 QAOA scheme for 3 qubits using p= 2.................. 16 4.1 One-dimensional beam target steering task at the CERN TT24-T4 beamline.................................. 19 4.2 Optimal QBM agent evaluation for the CERN TT24-T4 beam line. . 20 5.1 UMA class diagram for qmb-rl-steering . . . . . . . . . . . . . . . . . 21 5.2 D-Wave 3x3 Chimera graph, denoted C3. Qubits are arranged in 9 unitcells. ................................. 22 6.1 QAOA with p = 1, 2 compared to SQA. . . . . . . . . . . . . . . . . 27 6.2 RQAOA compared to SQA. . . . . . . . . . . . . . . . . . . . . . . . 28 6.3 QAOA and RQAOA compared to classical benchmark. . . . . . . . . 29 v CHAPTER 1. CLASSICAL REINFORCEMENT LEARNING 6 conventional Q-learning algorithm can be enhanced through the integration of Deep Learning, as we will see in the next section, with Deep Q-learning. 1.3 Deep Q-learning Deep Q-Learning (DQN)16 is a variant of Q-learning that uses deep neural networks to approximate the Q-function. The neural network is trained using a variant of supervised learning, where the input to the network is the current state of the environment and the outputs are the estimated Q-values for each possible action. The network is trained by minimizing the difference between the estimated Q-values and the target Q-values, which are updated during the interaction with the environment. DQN has been applied successfully to a variety of challenging problems in diverse areas and it is used as a classical benchmark (and inspiration) for FERL. Here we have the pseudocode for the full algorithm: Algorithm 2 Deep Q-learning with Experience Replay Initialize memory Dwith capacity N Initialize action-value function Qwith random weights θ Initialize target action-value function ˆ Qwith weights θ−=θ for episode 1, M do Initialize sequence s1={x1}and preprocessed sequence ϕ1=ϕ(s1) for t= 1, T do With probability ϵselect a random action atotherwise select at= argmaxaQ(ϕ(st), a;ϕ) Execute action atin emulator and observe reward rtand image xt+1 Set st+1 =st, at, xt+1 and preprocess ϕt+1 =ϕ(st+1) Store transition (ϕt, at, rt, ϕt+1) in D Sample random minibatch of transitions (ϕj, aj, rj, ϕj+1) in D Set yi=   rjif episode terminates at step j+ 1 rj+γmaxa′ˆ Q(ϕj+1, a′;θ−) otherwise Perform a gradient descent step on (yj−Q(ϕj, aj;ϕ))2with respect to the network parameters θ Every Csteps reset ˆ Q=Q end for end for We observe that in the last step ˆ Qis updated with Qevery Csteps. This is known as hard update. There is a variation known as a soft update (or Polyak update17) which consists of: ˆ Q=τQ + (1 −τ)ˆ Q(1.3.1) depending on τ, the soft update factor. Instead of having a huge update every C steps, we progressively update the function, improving the learning stability. Chapter 2 Quantum computing Quantum computing is one of the most promising and exciting emerging technologies that exist nowadays. Shor’s algorithm18 for integer factorization and Grover’s search algorithm19, as examples, describe powerful applications of these devices. The practical utility of these known applications is still out of our reach: we need much more powerful machines, the so-called fault-tolerant20 devices, which might be available in a few decades thanks to the fast development of Quantum Error Correction21. Current devices cannot still reach the high expectations that quantum computing promises. We currently find ourselves in the Noisy Intermediate-Scale Quantum (NISQ) era of quantum computers22; the devices available are quite noisy and circuits must have small depth so as to maintain the coherence of the quantum states. However, it is crucial for scientists and engineers working in this field to show the maximum potential of the devices that are available as of today. These demonstrations will lead to further investment for the development of these new technologies. One of the main goals in quantum computing is to demonstrate quantum advantage, which is proving that quantum computers can be much more efficient than conventional computers for some computational problems. The goal now is to start demonstrating quantum advantage with NISQ algorithms23. 2.1 Gate-based model fundamentals The basic unit of information in classical computers is the bit, which can only take values 0 or 1. In quantum information, the basic unit of information is the qubit, which can be defined as a superposition of the states 0 and 1, in the computational basis, with complex amplitudes: |q⟩=α|0⟩+β|1⟩=α βα2+β2= 1 α, β ∈C(2.1.1) 7 CHAPTER 2. QUANTUM COMPUTING 8 The squared modulus of this amplitudes represents the probability of collapsing the wave function onto either of the two states of the basis when measuring and have to fulfill the constraint of adding up to 1. Instead of using the complex amplitudes one can also use real angles θand ϕ: |q⟩= cos θ 2|0⟩+ sin θ 2eiϕ|1⟩θ∈[0, π], ϕ ∈[0,2π](2.1.2) This definition allows pure single qubit states to be represented on the surface of the Bloch sphere, a three-dimensional space where orthogonal qubit states point to opposite directions: Figure 2.1: Bloch sphere; a geometrical representation of a two-level quantum system. By Glosser.ca - Own work, CC BY-SA 3.0 In a quantum computer qubits are manipulated with quantum gates. These gates can be defined as operators acting on the Hilbert space H, describing the time evolution of a closed quantum system depending on the Hamiltionian H as U=e−iHt. These operators are reversible, and they can be represented as unitary matrices, i.e., for an operator Uit must hold that U†U=I. They can be applied to one qubit (single-qubit gates) or more (gates involving two or more qubits) at the same time. The most important single-qubit quantum gates are the ones defined by the Pauli matrices: σx=0 1 1 0, σy=0−i −i0, σz=1 0 0−1(2.1.3) The set formed by the Pauli matrices and the Identity matrix can be used as a basis to represent any single-qubit gate as a linear expression. The Hadamard gate can be written as H = (σx+σz)/√2: |0⟩H1 √2(|0⟩+|1⟩) Figure 2.2: Hadamard gate. H=1 √21 1 1−1(2.1.4) CHAPTER 2. QUANTUM COMPUTING 9 An important set of gates are the phase-shift gates, which generate a phase-shift between the two components of the qubit. The general formula for the phase-shift is: P(φ) = 1 0 0eiφ(2.1.5) Two known useful phase-shift gates are S and T, defined by: S=Pπ 2=1 0 0eiπ 2=1 0 0i=√σz(2.1.6) T=Pπ 4=1 0 0eiπ 4=√S=4 √σz(2.1.7) For universal quantum computation with pure states, entangling gates are necessary. The CNOT is the controlled version of the Pauli-X gate (also know as NOT gate), and it’s a maximally entangling gate. One qubit acts as a control and the other as a target. |a⟩•|a⟩ |b⟩ |a⊕b⟩ Figure 2.3: CNOT gate. CNOT =    1 0 0 0 0 1 0 0 0 0 0 1 0 0 1 0     (2.1.8) It is proven that gates H, S and CNOT (known as the Clifford set), alongside the T gate form a fault-tolerant universal quantum basis24, i.e. approximates arbitrarily well any given unitary. Quantum algorithms are typically represented using quantum circuits, which are analogous to classical electronic circuits. In these circuits, each qubit is depicted as a line. As an illustrative example, Figure 2.4 depicts the quantum teleportation protocol25. This particular circuit involves three qubits, with the objective of teleporting the first qubit to the third qubit while using the second qubit solely for performing operations (referred to as an ancilla qubit). The circuit incorporates two Hadamard gates and two CNOT gates. To successfully execute the teleportation process, it is essential to measure the first two qubits, which subsequently leads to their destruction, conditioning the application of the Pauli-X and Pauli-Z gates on the teleported qubit. CHAPTER 2. QUANTUM COMPUTING 10 |ψ⟩•H• |0⟩H• |0⟩X Z |ψ⟩ Figure 2.4: Quantum teleportation circuit. 2.2 Adiabatic quantum computing Adiabatic quantum computing (AQC)26 is a specialized approach to quantum computing that leverages the concept of quantum annealing. In quantum annealing, a quantum system is slowly evolved from an initial Hamiltonian that represents a simple problem to a final Hamiltonian that encodes the desired problem. The hope is that the system remains in its ground state throughout this adiabatic evolution, thus providing a solution to the problem encoded in the final Hamiltonian. H(t) = A(t)Hi+B(t)Hf(2.2.1) where A(t), B(t) : [0, T]→Rare such that A(0) = B(T) = 1 and A(T) = B(0) = 0. Provided that the evolution is slow enough, the adiabatic theorem ensures that, starting from the ground state of H(0) = Hi, the system will end up in the ground state of H(T) = Hf. Quantum annealing (QA) is often compared to simulated annealing, a classical algorithm used for optimization. In simulated annealing, the exploration of the solution space is accomplished through iterative transitions between states, where the probability of these transitions is determined by a specific parameter called the temperature. While both approaches involve annealing processes, quantum annealing exploits quantum effects such as superposition and entanglement to explore multiple states simultaneously. This quantum advantage can potentially lead to more efficient optimization in multiple cases. The computational simulation of the quantum annealing process, which will be performed in this project, is known as simulated quantum annealing (SQA), and it has to be discerned from simulated annealing. It is important to note that adiabatic quantum computing is not capable of universal quantum computation. Unlike gate-based quantum computing architectures, which can perform any computation through a sequence of quantum gates, AQC is specifically designed for solving optimization problems. While this limitation restricts the scope of problems that can be tackled using AQC, it also allows for a simplified and more efficient implementation. The primary application of adiabatic quantum computing is in solving optimization problems. By mapping the problem to be solved onto the final Hamiltonian of the quantum system, AQC aims to find the ground state that corresponds to the optimal solution. This approach has been shown to be particularly effective for problems that can be formulated as finding the minimum energy configuration in a CHAPTER 2. QUANTUM COMPUTING 11 complex energy landscape. Comparisons between quantum annealing and Quantum Approximate Optimization Algorithm (QAOA) are also relevant in the context of adiabatic quantum computing27. QAOA is a variational algorithm that uses a parameterized quantum circuit to approximate the ground state of a given problem. While both quantum annealing and QAOA aim to find optimal solutions, they apparently employ different approaches. Quantum annealing leverages adiabatic evolution, while QAOA utilizes quantum gates and measurements. The choice between these two approaches depends on the specific problem at hand, as different problems may benefit from one method over the other. Even though QAOA has been defined as a Trotterized or truncated approach of QA5, there are results showing that QAOA is able to deterministically find the solution of specially constructed optimization problems in cases where QA fails27. CHAPTER 2. QUANTUM COMPUTING 12 Chapter 3 Quantum Machine Learning Quantum Machine Learning (QML)28 combines the principles of quantum computing and machine learning to perform artificial intelligence tasks such as regression, classification and clustering. Unlike classical machine learning algorithms, QML allows for the exploitation of quantum mechanical phenomena, such as superposition and entanglement, to process data in ways that are not possible with classical computing. This leads to a potential speedup in processing time for some tasks and the ability to tackle new types of problems in areas such quantum chemistry and finance. optimization problems QML algorithms can be classified according to Figure 3.1 by the type of data it is using (quantum or classical) and whether it is performed in a conventional computer or on a quantum computer. When an algorithm uses both classical and quantum resources it is known as an hybrid algorithm (e.g. QAOA). Figure 3.1: Quantum machine learning diagram. By Maria Schuld - Own work, CC BY-SA 4.0 3.1 Quantum Boltzmann Machine A Boltzmann machine (BM) is a type of generative stochastic artificial neural network. It consists of a set of visible and hidden (latent) variables. They are denoted by vi, with i∈V, and by hj, with j∈H, respectively. Visible variables typi13 CHAPTER 3. QUANTUM MACHINE LEARNING 14 cally serve as the inputs and outputs, and hidden variables are used to increase the expressivity of the model. For the FERL algorithms discussed in this paper, the BM’s topology is chosen to be clamped, building upon the research in Levit et al.3. In a clamped topology the visible variables are not part of the network; they enter the system’s energy function as biases or self-couplings to the specific hidden variables assigned to them. These bias terms are used in a weighted sum, where the weights are given by the coupling strengths between visible and hidden variables. Bias terms are linear in the state of the hidden variable with which they are associated and couplings between hidden variables, such as wjk ∈R,j, k ∈H, correspond to quadratic contributors to the energy function. For reinforcement learning, the visible variables are composed of the vectors of the state-action pair v= (s, a)∈RdS+dA, where dSand dAare the dimensions of the state and action space. Figure 3.2: Clamped Boltzmann Machine used in FERL. By Schenk et al. A Quantum Boltzmann Machine (QBM) is a quantum computing algorithm that is inspired by the classical Boltzmann machine, which was first proposed by Amin et al. in 201629. It can be represented by a physical system of coupled qubits in the presence of a purely transverse magnetic field, here pointing along the x-axis. Each node of the BM can assume a spin up or down state following a certain probability distribution. The system’s energy states are described by the Hamiltonian of the transverse-field Ising model30: H(v) = −X i∈V, j∈H ωijviσz hj−X j,k∈H ωjkσz hjσz hk−ΓX j∈H σx hj(3.1.1) where σx hj and σz hjare the Pauli spin matrices (introduced in Chapter 2) acting on the hidden node hjfor the xand zdirections, respectively, and Γdenotes the transverse magnetic field strength. The transverse field introduces quantum fluctuations which enable, among others, quantum tunnelling [14]. When spin states are measured along the zcoordinate, information about the spin xcomponents becomes lost. This limitation can be addressed through the application of replica stacking. This approach was introduced by Levit et al.3to extend the transverse-field Ising model using the Suzuki-Trotter expansion. By implementing this method, it becomes possible to depict the non-zero transverse- CHAPTER 3. QUANTUM MACHINE LEARNING 15 field Ising model as a classical Ising model in a dimension that is one higher. The effective Hamiltonian for this extended representation is as follows: Heff(v) = −1 NrX l=1    X j,k∈H ωjkhj,lhk,l +X i∈V, j∈H ωijvihj,l   −ω+ X j∈H Nr X l=1 hj,lhj,l+1 +X j∈H hj,1hj,Nr! (3.1.2) with ω+=1 2βlog hcothΓβ Nri, where Nris the number of replicas, βthe inverse temperature, and hj,l denotes the hidden node with index jin replica l. 3.2 Quantum Approximate Optimization Algorithm (QAOA) The Quantum Approximate Optimization Algorithm5(QAOA) is a quantum machine learning algorithm designed for gate-based model computers to produce approximate solutions for combinatorial optimization problems. The algorithm depends on an integer p≥1and the quality of the approximation improves as pis increased. The quantum circuit that implements the algorithm consists of unitary gates whose locality is at most the locality of the objective function whose optimum is sought. The depth of the circuit grows linearly with ptimes (at worst) the number of constraints. First, we have to define a cost Hamiltonian HC, which depends on the optimization problem we aim to solve, and a mixer Hamiltonian HM, which is a key element for the algorithm to evolve successfully. Once these two Hamiltonians are defined, we construct the transformations e−iγHCand e−iαHM. The parameters αand γare the parameters being optimized in the algorithm through classical methods. The set of both parametrized unitary transformations defined by the cost and the mixer Hamiltonians form a layer of the QAOA algorithm (highlighted in Figure 3.3). The number of layers is determined by the parameter pintroduced earlier, which is an hyperparameter. The final circuit can be defined as: U(γ, α) = e−iαpHMe−iγpHC... e−iα2HMe−iγ2HCe−iα1HMe−iγ1HC(3.2.1) Once the algorithm is defined, we can make the observation that the mixer operator is necessary to prevent the algorithm from getting stuck in one of the HC eigenstates, mixing the probabilities of all possible solutions. For that purpose, the mixer Hamiltonian has to be defined that it does not commute with the cost Hamiltonian (HCHM=HMHC). The original QAOA uses the X-mixer (applying the parametrized σXto each qubit) but there are others, and it is plausible to evaluate which mixer works best depending on the problem in question31. CHAPTER 5. IMPLEMENTATION 22 implemented using the Sqaod43 library, and the QAOA solver, which is explained in the following section. The Boltzmann machine couplings are following the Chimera graph topology of the D-Wave 2000Q QPU that can be observed in Figure 5.2. This architecture is comprised by sets of connected unit cells. In Figure 5.2, the unit cell is rendered as a column, but it could be also rendered equivalently as a cross. Each unit cell has eight qubits (n_nodes_per_unit_cell=8): four horizontal qubits connected to four vertical qubits via couplers (bipartite connectivity), e.g. observing at the first unit cell, each qubit from 0 to 3 in the first column is connected to every qubit 4 to 7 in the second column. Unit cells are ordered vertically (rows represented by n_rows) and horizontally (columns represented by n_columns) with equivalent qubits connected, creating a lattice of sparsely connected qubits, e.g. observing at the central unit cell, qubits 32 to 35 are each connected to the nodes that occupy its equivalent relative position in the unit cells in the same column and neighbouring rows (nodes 8 to 11 and 56 to 59) and, likewise, qubits 36 to 39 are each connected to its analogous node in the unit cells in the same row and neighbouring columns (nodes 28 to 31 and 44 to 47). Taking all the connections into account, Chimera qubits are considered to have a nominal length of 4 (each qubit is connected to 4 orthogonal qubits through internal couplers) and degree of 6 (each qubit is coupled to 6 different qubits). Figure 5.2: D-Wave 3x3 Chimera graph, denoted C3. Qubits are arranged in 9 unit cells. By D-Wave In our algorithm we choose n_columns=2 and n_rows=1, obtaining a n_unit_cells=2 and multiplying by the number of nodes per cell we get a total of n_nodes=16. Observing Figure 5.2, our structure could be represented by the nodes comprised between 0 and 15. Our algorithm is conditioned by a significant amount of hyperparameters (most of them introduced in previous chapters), that need to be determined before the training, and can be later tuned in order to increase its performance. Our aim is to find and work with the best parameters. These parameters can be divided CHAPTER 5. IMPLEMENTATION 23 into three main categories: exploration parameters, learning parameters and solver parameters. We start with the first two categories, indicating the parameters used in our implementation, which were previously tuned using Optuna44: •Exploration parameters. They are the parameters directly related to the reinforcement learning process, controlling how the agent interacts with the environment (see Chapter 1): -exploration_epsilon: it is a tuple of initial and final ϵ-greedy factor introduced in Section 1.1, decaying linearly over the exploration_fraction. We fix it to (1.0, 0.), meaning we start with a random policy, and we finish with a completely greedy policy. -exploration_fraction: the fraction of training period with ϵ-greedy policy RL parameter. We fix it to 0.766. -max_steps_per_episode: maximum allowed iterations for the agent to interact with the environment per episode. We fix it to 15. •Learning parameters. They are the parameters related to the learning process, controlling the way the agent learns through its interactions with the environment: -learning_rate: it is a tuple with the learning rate at the start and the end, with a linear decay, to update the coupling weights of the Chimera graph. It is fixed at (0.04778, 0.000201). -small_gamma is the reward discount factor γfor cumulative future rewards introduced in Equation 1.2.2. It is fixed at 0.756. -replay_batch_size: the batch size i.e. the number of experiences that we base our training on in every training iteration, sampled from replay buffer. It is fixed at 32. -target_update_frequency: determines the number of iterations after which we update the target Q-function following soft-update rule. It is set to 1 (every iteration). -soft_update_factor: the τparameter for the Polyak update of the target Q-function net weights introduced in Equation 1.3.1. It is set at 0.4008. The last category of parameters depend on the solver we are using. We will explain the parameters regarding the SQA (and equivalently QPU) solver, which will not be present when using the QAOA solver. •Solver parameters. For SQA, the parameters determine the evolution of the simulated annealing process: -n_replicas: number of replicas in the 3D extension of the Ising spin model (Trotter slices) to account for non-zero transverse field Γat the end of the annealing process (see Figure 1 in Levit et al.3) and was introduced in Equation 3.1.2 as Nr. We set it to 1. CHAPTER 5. IMPLEMENTATION 24 -n_meas_for_average: number of times we run an ’independent spin configuration sampling process, which will be used for the SQA and eventually to calculate the average effective Hamiltonian of the system. We set it to 16. -n_annealing_steps: number of steps that one annealing process should take. We fix it to 180. -big_gamma: a tuple indicating strength of the transverse field Γintroduced in Equation 3.1.1 decaying from the first to the second value over the course of the quantum annealing process. We set it to (41.41, 0.). -beta: the inverse temperature introduced in Equation 3.1.2 (note that this parameter is kept constant in SQA unlike simulated annealing). We fix it to 0.746. In this list of fixed parameters we have not included the number of training steps (n_steps), which belongs to the learning parameters category. In our executions, we will scan our algorithm versions through an array of values. This scan will help us compare different solvers and its different variations, to compare the performance in terms of training steps, and hence find the algorithm which can perform better using less steps. 5.2 QAOA implementation as solver To implement the QAOA algorithm described in Section 3.2 we use IBM’s Qiskit45 library, which has a function that directly implements the algorithm and allows us to change the configuration depending on the inputs. 1class QAOA: 2def __init__ ( self , n_nodes : int = 16 , solver : str =’QAOA ’ 3p: int = 1) -> None : 4 5self. n_nodes = n_nodes 6self.p = p 7 8backend = qsk .Aer . get_backend (’aer_simulator ’) 9backend . set_options ( device = ’CPU ’) 10 sampler = BackendSampler ( backend ) 11 optimizer = COBYLA () 12 qaoa_problem = qsk . algorithms . minimum_eigensolvers . QAOA ( sampler , optimizer , reps = self.p) 13 14 self .solver = MinimumEigenOptimizer ( qaoa_problem ) 15 16 def _reformulate_qubo (self , qubo_dict ) -> QuadraticProgram : 17 qubo_problem = QuadraticProgram () 18 19 # Define all the " qubits " / binary variables 20 for iin range ( self . n_nodes ): 21 qubo_problem . binary_var (’x’ +str(i)) 22 23 linear_terms = np . zeros ( self . n_nodes ) CHAPTER 5. IMPLEMENTATION 25 24 quadratic_terms = {} 25 for (j, k), w in qubo_dict . items () : 26 if j == k: 27 linear_terms[j] = w 28 else: 29 quadratic_terms[(’x’ +str (j), ’x’ +str(k))] = w 30 31 qubo_problem . minimize ( linear = linear_terms , quadratic = quadratic_terms ) 32 33 return qubo_problem 34 35 def sample (self , qubo_dict : Dict , n_meas_for_average : int, 36 *args , ** kwargs ) -> np . ndarray : 37 38 qubo_problem = self . _reformulate_qubo ( qubo_dict ) 39 spin_configurations = [list( self . solver . solve ( qubo_problem ) .x)] 40 41 # Convert to np array and flip all the 0s to -1s 42 spin_configurations = np.array(spin_configurations) 43 spin_configurations[spin_configurations == 0] = -1 44 45 spin_configurations = spin_configurations.reshape( 46 (1 , 1, self . n_nodes )) 47 48 return spin_configurations Listing 5.1: QAOA solver implementation using Qiskit The algorithm is implemented with the minimum_eigensolvers.QAOA() function, which is given three input arguments: sampler,optimizer and reps. For the sampler we have used the aer_simulator backend, whose default behavior is to mimic the execution of an actual device. The classical optimizer chosen is the COBYLA. Finally the parameter reps represents the parameter pintroduced in Section 3.2, which is the number of times we repeat the algorithm, adding more layers of the cost and mixer Hamiltonians, which will be our main hyperparameter. Once the qaoa_problem is set we have to create the solver using the MinimumEigenOptimizer() function from the qiskit-optimization package. In our original implementation, the QBM is represented as a Quadratic Unconstrained Binary Optimization (QUBO) problem, since this formulation can be directly applied to quantum annealers. Finding the solution to a QUBO is equivalent to finding the ground state of a corresponding Ising Hamiltonian, as the one shown in Equation 3.1.1 to define the QBM46. Ising formulation uses Ising space si∈ {−1,1} and QUBO uses Boolean space xi∈ {0,1}, and they can be transformed into one another through a linear transformation: xi=1 + si 2 The reformulate_qubo function converts the input QUBO dictionary into the QuadraticProgram() Qiskit class. After that, using the solve method in our CHAPTER 5. IMPLEMENTATION 26 solver, it first wraps the translation to an Ising Hamiltonian again, which is the formulation used for the QAOA. After that, it calls the MinimumEigensolver and obtains the translation of the results back to an OptimizationResult, whose optimal value can be obtained using .x. We have now completed the evaluation of the QBM, and its parameters can be optimized through the rule defined in Equation 3.3.2. To now implement RQAOA, also introduced in the past Section 3.2, we have to introduce some modifications to our code. 1class QAOA: 2def __init__ ( self , n_nodes : int = 16 , solver : str =’RQAOA ’, 3p: int = 1, min_num_vars: int = 100) -> None: 4 5self. n_nodes = n_nodes 6self.p = p 7self.min_num_vars = min_num_vars 8 9backend = qsk .Aer . get_backend (’aer_simulator ’) 10 backend . set_options ( device = ’CPU ’) 11 sampler = BackendSampler ( backend ) 12 optimizer = COBYLA () 13 qaoa_problem = qsk . algorithms . minimum_eigensolvers . QAOA ( sampler , optimizer , reps = self.p) 14 15 self . solver_internal = MinimumEigenOptimizer ( qaoa_problem ) 16 self . solver = RecursiveMinimumEigenOptimizer ( self . solver_internal , min_num_vars = self . min_num_vars ) 17 18 (...) Listing 5.2: RQAOA solver implementation using Qiskit The RecursiveMinimumEigenOptimizer() applies a recursive optimization on top of OptimizationAlgorithm. This optimizer can use MinimumEigenOptimizer as an optimizer that is called at each iteration. The algorithm is introduced in Bravyi et al.47. Chapter 6 Results In this section we present the results that we have obtained following the procedure described in the previous Chapter 5 for the FERL algorithm using QAOA and SQA solvers applied to the 1D proton beam line environment at CERN described in Chapter 4. Our first results are presented in Figure 6.1. Each of the curves has been obtained by averaging 10 repetitions. It represents the optimality of the policy for increasing number of training steps. The optimality is calculated comparing the policy obtained using the algorithm to the optimal algorithm, which is known for this environment, using a total of 200 sample states. Figure 6.1: QAOA with p = 1, 2 compared to SQA. 27 CHAPTER 6. RESULTS 28 In this figure we present the results obtained using QAOA for p= 1,2compared to SQA. We observe that the best results are obtained with the SQA, reaching higher optimality values for lower number of steps, and reaching 100% for n_steps = 30. Both QAOA versions are similar in performance in the first stages, but for n_steps = 30, for p= 1 it approaches full optimality (reaching 100% at n_steps = 50) while for p= 2 it falls to approximately 80%. We can affirm that for this specific problem QAOA with p= 1 performs better than with p= 2. This affirmation might not seem in line with what we stated in Chapter 3: the quality of the approximation improves as pis increased. However, when increasing pwe are growing the complexity of the problem to solve, increasing the number of local minima and difficulting the optimization48. For problems with higher complexity and higher number of qubits, a higher value of pmight be needed, but for this case p= 1 is enough. During the run, we can qualitatively observe that the run-time per execution increases with every time-step until it reaches the step 32, which is the replay_batch_size (see Section 5), where it stabilizes. We also observe that the computational cost and thus the execution time for p= 2 is considerably higher than for p= 1. Figure 6.2: RQAOA compared to SQA. Next, in Figure 6.2 we compare the RQAOA solver with the SQA. For this case, for each solver we have averaged a total of 50 runs. We observe that the performance for both solvers is similar, though we can argue that RQAOA outperforms SQA, since it reaches 100% optimality before at n_steps = 20, before than SQA. CHAPTER 6. RESULTS 29 Figure 6.3: QAOA and RQAOA compared to classical benchmark. Finally, in Figure 6.3, we combine the information obtained in the two previous for the QAOA with p= 1,2comparing it to the classical solver benchmark provided by Qiskit. We clearly observe how RQAOA outperforms original QAOA, and that it adjusts well to the classical benchmark, which serves as an upper limit. It has already been demonstrated that RQAOA outperforms the original QAOA for other problems such as the MaxCut49. CHAPTER 6. RESULTS 30 Chapter 7 Conclusion Throuhout this project we have been able to successfully describe the FERL algorithm giving a significant theoretical background in reinforcement learning and quantum computing to understand its main components, applying it to the 1D proton beam line at environment at CERN. We have proven that the QAOA algorithm can be effectively implemented as a solver for the QBM in FERL. This implies that the algorithm can be adapted to a gate-based universal quantum architecture instead of relying on a quantum annealer. We have observed that for this specific environment, given its low complexity, it is better to keep the repetition parameter pat its minimum value of 1. We have also observed that RQAOA outperforms the original QAOA and that, at this stage of simulation, it can outperform the SQA solver. Further work on this implementation would include an extension of the simulations including experimental noise and, eventually, the implementation of the algorithm into real quantum hardware. Additionally, another goal is to evaluate this algorithm with a broader perspective, instead of directly focusing on the beam-steering task, trying to make these QRL approaches work in prototypical RL benchmark problems, such as the cart-pole problem. The QAOA solver can now be seen as a building block for further algorithms, e.g. the hybrid actor-critic scheme in development at CERN, whose performance could be increased by the use of this new solver and compared again to the classical DDPG, allowing to efficiently apply QRL to control problems with a continuous action space. 31