Optimal quantum control via genetic algorithms for quantum state engineering in driven-resonator mediated networks

We employ a machine learning-enabled approach to quantum state engineering based on evolutionary algorithms. In particular, we focus on superconducting platforms and consider a network of qubits—encoded in the states of artificial atoms with no direct coupling—interacting via a common single-mode driven microwave resonator. The qubit-resonator couplings are assumed to be in the resonant regime and tunable in time. A genetic algorithm is used in order to find the functional time-dependence of the couplings that optimise the fidelity between the evolved state and a variety of targets, including three-qubit GHZ and Dicke states and four-qubit graph states. We observe high quantum fidelities (above 0.96 in the worst case setting of a system of effective dimension 96), fast preparation times, and resilience to noise, despite the algorithm being trained in the ideal noise-free setting. These results show that the genetic algorithms represent an effective approach to control quantum systems of large dimensions.


I. INTRODUCTION
Quantum state engineering is an essential enabling step for a variety of quantum information tasks, including the initialization of quantum simulators [1], the loading of classical data for quantum-enhanced analysis [2], or the generation of resourceful states in quantum communication networks [3].In particular, quantum entangled states, which embody a stark departure from classicality, often provide the main resources towards quantum advantage [4].Therefore the preparation of entangled states in multi-node quantum networks presents a key challenge for realising quantum protocols on near term devices.
A successful approach towards state and resource generation consists of two steps: (i) engineer a suitable timedependent Hamiltonian with tunable parameters, as allowed by current experimental capabilities in a given physical platform; (ii) find and implement proper temporal dependencies (pulse shapes) for these parameters by invoking optimal quantum control techniques [5,6].This results in tailor-made control schemes pertinent to the specific platform at hand.
One platform at the forefront of engineering flexible multinode quantum networks, to which such approach has been successfully applied, is that of superconducting quantum circuits [7][8][9][10].So far superconducting architectures have been employed to realise two-qubit gates using frequencytunable [11][12][13][14] and microwave-driven [15,16] artificial atoms.In addition, demonstration of the coupling between artificial atoms and microwave resonators [17] opened the door for resonator mediated two-qubit gates [18][19][20][21] and provided an alternative platform to study cavity quantum electrodynamics [22] leading to the field of Circuit QED [23].Whilst extremely flexible in their design, it has been shown however that operating these superconducting systems with a reduced level control is not only desirable, but necessary in some cases [24,25].Thus finding optimal control protocols that utilize a limited but effective level of control is of practical interest.
Optimal control of quantum systems has yielded a range of new methods inspired, in part, by the development of modern machine learning methods.Specifically, neural-networkbased reinforcement learning methods [26,27] have been recognised as useful tools to study quantum systems [28] in a variety of contexts including state transfer [29][30][31][32], quantum thermodynamics [31], circuit architecture search [33], control of dissipative systems [34] and Adaptive Quantum Enhanced Metrology [35].Reinforcement learning techniques have proven particularly suitable for control problems of increasing dimension when compared to more standard techniques [36].However, for the most part these techniques cannot be used as closed-loop optimization schemes and therefore are of limited use for optimization on physical quantum systems [37] and have relatively poor convergence guarantees which necessitate, often expensive, hyper parameter tuning steps.This motivates one to consider whether these control problems can be tackled using comparatively less sophisticated techniques that scale well but have greater ease of implementation.With this in mind we consider evolutionary strategies, which have already been proposed as a scalable substitute to reinforcement learning methods [38] and used as an alternative to gradient based parameter updates in both deep reinforcement learning [39] and quantum reservoir computing [40].In fact, Natural Evolution Strategies (NES) have already been proposed in the context of quantum control [35,41] resulting in the the realisation of fast, high-fidelity single-shot three-qubit gates in frequency tunable superconducting systems [42,43], while extension to four qubit systems has also been tackled effectively using more standard global optimization techniques [44].More recently these NES techniques have been utilized to address optimal annealling schedules in spin systems [45] and optimal transport of Majorana fermions in superconducting nanowires [46].The aim of this work is to investigate the potential of these techniques in the specific context of direct state preparation in a drivenresonator mediated, reduced-control multi-qubit network.
Here we consider a register of qubits coupled via a common resonator and operated in a regime which facilitates a re-duced level of control, and employ a genetic algorithm to find optimal pulse sequences to drive their dynamics.In order to illustrate our approach, we present efficient control schemes for preparing entangled three-and four-qubit states, including GHZ, Dicke, and graph states, and assess performance against relevant decoherence sources finding the thresholds that limit the quality of our results.The generation of high-quality states, in short times compared to typical multi-qubit circuit timescales, is thus demonstrated while identifying fully the sequence of driving pulses to use.Our approach contributes to the growing argument that a hybrid take to quantum control -that mixes machine learning and optimal control -is a viable route to the engineering of crucial resources for quantum information processing.
The remainder of the manuscript is organised as follows.In Sec.II we present the specifics of the system considered and formalise the Hamiltonian.In Sec.III, we present an overview of the Continuous Genetic Algorithm employed then follow by formalising the algorithm for the problem of quantum control in section IV.Sec.V focuses on the presentation of the resulting control schemes for the preparation of three-qubit GHZ and Dicke states as well as a specific instance of a fourqubit graph states.In Sec.V C, the effect of decoherence is investigated.Sec.VI offers our closing remarks and a forward look.

II. SYSTEM
We consider a system composed of N identical and noninteracting qubits coupled via a common single-mode driven resonator.We model such system with the Hamiltonian H = H 0 + H int + H d , where In writing these Hamiltonians, we have assumed the rotating wave approximation and used units such that = 1.Here a denotes the annihilation operator of the resonator mode, whereas σ + j = |1 0| and σ − j = |0 1| are the raising and lowering operators of the j-th qubit.Moreover, ω c denotes the resonator frequency, ω j the transition frequency of the j-th qubit and g j its coupling strength with the resonator.The driving amplitude (assumed to be real) and the carrier frequency of the drive are indicated with ξ and ω d respectively.
The use of a resonator to mediate two-qubit gate interactions is well studied.For example, in the context of superconducting systems the Resonator-Induced Phase (RIP) gate [18][19][20] utilizes two qubits dispersively coupled to a common microwave resonator.In such a dispersive regime, the coherent qubit-resonator interaction term becomes negligible leading to an effective qubit-qubit interaction term.Whilst this can be beneficial for protection against decoherence there are two main drawbacks: 1) The dispersive coupling precludes real photon processes between the resonator and each qubit which necessitates the use of local qubit drives in order to introduce energy into the system.This subsequently makes this regime adverse to scaling to larger numbers of qubits.2) The large detunings involved result in small effective qubit-qubit coupling strengths at the expense of gate operation time.On the other end of the spectrum, resonant regimes have been considered to selectively tune qubits in and out of coupling with a common resonator [47].In addition, external microwave drives are also commonly used such as with the the cross-resonance gate [15,16].In contrast to such previous works, we assume in H a fully resonant regime which allows us to adapt a limited control scheme, in that the local qubit drives can be replaced by a single resonator driving term where the amplitude of this drive and the qubit-resonator couplings are tunable.This facilitates the possibility of multi-qubit interactions that operate on timescales that are much shorter than typical multi-qubit gate times.
Taking the interaction picture with respect to the frequency of the harmonic mode and using H = e iRt He −iRt − R, where where δ j = ω j −ω c and δ d = ω d −ω c are the detunings between the j th qubit and the harmonic mode, and between the drive frequency and the harmonic mode respectively.The utility of this transformation becomes apparent when one assumes full resonance between the qubit transition frequencies, drive frequency and resonator frequency -namely, ω c = ω d = ω j , ∀ j = 1, ..., N. The Hamiltonian thus takes the simpler form which comprises N + 1 terms embodying an equal number of controls, where we assume each of the qubit-cavity couplings g j (t) and the drive amplitude ξ(t) to be time-dependent controllable parameters.
The above Hamiltonian governs the dynamics of the system via the time-dependent Schrödinger equation.Thus if the system is initially in some state |ψ(t 0 ) , then the time-evolved state of the system at any future time t > t 0 is given by [48] The goal here is to find the optimal functional forms for g j (t) and ξ(t), such that the system is dynamically steered in some desired way.Specifically, we are interested in state preparation within the qubit subspace, so we first determine the reduced state of the qubit network and work to find g j (t) and ξ(t) such that the fidelity is maximized, where |ψ σ is the target state of interest.It is worth stressing that our approach would work equallymutatis mutandis -with mixed target states.

III. THE CONTINUOUS GENETIC ALGORITHM
Evolutionary strategies are a class of direct search optimization techniques, drawing inspiration from Darwinian evolution, that have recently been proposed as a viable substitute for gradient based parameter optimization in Neural Networks [39] and quantum reservoirs [40], as well as a scalable alternative to Reinforcement learning techniques [38].Of particular interest to continuous-control problems is the so-called "Continuous Genetic Algorithm" (CGA) [49], which generalises the more traditional Discrete Genetic algorithm to allow for continuous parameter values.
Consider an optimization problem with N var parameter variables p i .We call a specific instance of these parameters, C = [p 1 , p 2 , ..., p N var ], a Chromosome.This embodies one proposed solution to the optimization problem.One then defines a Fitness function, f (C) = f (p 1 , p 2 , ..., p N var ): R N var → R.This function will be determined by the optimization task under consideration and will assign a numerical score to each proposed solution (chromosome) depending on its usefulness (fitness) with respect to the specified task at hand.In the name of generality we assert that p i ∈ [−1, 1], where the parameter values are suitably scaled within the fitness function.We call each new iteration of the algorithm a Generation where the algorithm proposes several chromosomes to make up a Population.The algorithm can be qualitatively summarised using the following steps: 1. Initialization.Define the fitness function (determined by the optimization task) and fix the hyper-parameters for the genetic algorithm.One then generates an initial population with N pop chromosomes, typically randomly, which acts as the zero-th generation of the algorithm.
2. Natural Selection.The fitness of each chromosome is assessed with a call to the fitness function.The population is then ranked based on these fitness values and the N survive highest scoring chromosomes are chosen to survive, while the rest are discarded to be replaced by new offspring chromosomes in the next generation.

3.
Pairing.From the N survive surviving chromosomes chose N parents pairs of parent chromosomes in order to produce (N pop − N survive ) offspring to repopulate.Parent chromosomes are sampled probabilistically based on their relative fitness, such that the fittest chromosomes reproduce more frequently, where cloning (using the same chromosome for both parents) is disallowed.

4.
Mating.Here, each pair of parent chromosomes is combined in some manner to produce enough offspring to replenish the population.In the simplest case, one can consider cutting the parent chromosomes in half and producing two offspring, made up of the two possible combinations of these halves.Alternatively, one could pick indices along the chromosomes at random, and generate two offspring by first copying the parent chromosomes, then swapping the parameter values corresponding to these indices between the two copies.This direct transfer of parameter values to generate offspring is the simplest and most obvious way of mating two parent chromosomes.However this approach simply provides offspring chromosomes made up entirely of parameter values that were already present in the previous generation and, as such, introduces no new "genetic material" in the form of novel parameter values for a given parameter.To tackle this we can implement a combination of the form This allows us to introduce an element of exploration by allowing not only parameter values present in the previous generation but also any continuous value in between, modulated by random variable β ∈ [0, 1].

5.
Mutation.This final step introduces further exploration into the search by choosing, at random, a number of elements within each chromosome to be replaced with a new random value.The rate of this mutation is set by a parameter 0 < α < 1, which determines the number of indices to be targeted relative to N var .This mutation step is applied to the entire population except for the single chromosome with the highest fitness, which is known as Elitism.This is important to ensure theoretical guarantees of convergence.Specifically it ensures that the maximum achieved fitness is always at least maintained in new generations.
The algorithm repeats steps 2-5 until convergence or an acceptable level of fitness has been achieved.A key point to note here is that the CGA is a stochastic optimization process so with finite time the convergence to any absolute optimal strategy is never guaranteed, however on average one would expect an increase in performance as the computational time grows.The aforementioned hyper-parameters associated with CGA then are: population size N pop , number of survivors N survive to keep in each iteration, number of parental pairs to mate N parents and the mutation rate α.Also, as discussed above, we have some freedom in how we implement the mating procedure.

IV. CONTINUOUS GENETIC ALGORITHMS FOR OPTIMAL QUANTUM CONTROL
In order to apply the CGA we first must formulate the optimal control problem in a suitable manner.As is common practice, we discretize the evolution time into T time intervals of equal duration τ -which is manually chosen -and assume the functional form for each control to be defined by its values at the T + 1 times t = 0, τ, 2τ, ..., T τ (connecting the latter with a simple tanh function, similar to the use of piecewise error functions [42,43].See Appendix A for specifics).Therefore, given that each control function is completely defined by these T + 1 values, and there are N + 1 controls in total as per Eq. ( 5), then the total number of parameters that we have to optimize over is N var = (N + 1)(T + 1) (which as said is the length of the chromosomes).We then need to outline how we asses the fitness of each chromosome and in doing so define the optimization task.

A. The Fitness Function
As said, each chromosome is a sequence of N var parameter values, and there are N + 1 control pulses.Therefore we first break each chromosome up into N + 1 sequences where the first T + 1 parameters determine the first control pulse and so on, as in Fig. 1.Then the parameter values corresponding to each control are scaled accordingly where each of the qubitresonator couplings assume a common range g j ∈ [−g 0 , g 0 ], and the drive amplitude is in the range ξ ∈ [−ξ 0 , ξ 0 ].The scaled parameter values are then used to build the control pulses as in Appendix A, which fully define the time dependent Hamiltonian H(t) in Eq. 5.After defining an initial state |ψ(t = 0) the system evolves unitarily according to H(t) (as in Eq. 6), where we keep track of the state of the system at every intermediate time between t = 0 and t = τT .For each state in this history of states, we trace out the cavity system to obtain a state history for the qubit subspace, ρ Q (t).Next we calculate F (ρ Q (t), σ), which is the time dependent fidelity induced by the specific control pulses.In order to assign a numerical fitness value to the chromosome we simply take the maximum fidelity achieved in the qubit subspace throughout the induced dynamics, i.e In fact, the function actually used is a slight variation of that presented above necessitated by the specific details of the simulation and the ability to "extract" the state with the highest fidelity (see Appendix C for details).
V. RESULTS

A. Analysis of performance: case studies
Below we outline the optimal control schemes found when applying the CGA approach to prepare genuinely entangled 3 and 4 -qubit states from completely separable initial states.In doing so this acts as a proof of principle -for both the Hamiltonian H, and the optimization method -with respect to entanglement generation in general and state preparation specifically.Below we consider the physically reasonable maximum coupling and driving g 0 , ξ 0 = 2π 200MHz, and total time-scales τT ≤ 10ns (see appendix B for full details).An important point to make is that such timescales are considerably shorter than the typical gate times associated with two-qubit operations on superconducting platforms, which are commonly of the order of 100ns [8,9].This means that each of the following protocols operate on a much faster timescale than their circuit counterparts which prepare the same state.We specifically consider 3 states: GHZ state, three-qubit Dicke states with two excitation, and the 4 qubit "Box" Cluster state.
a. GHZ state.GHZ states are relevant genuinely tripartite entangled states [50] defined by We assume the system to initially be in the state |ψ 0 = |010 Q |0 cav .We use |010 as the initial state of the qubit network in this specific instance, as opposed to simply the global vacuum, since the latter has non-zero fidelity to the target GHZ state and thus starts the optimization in at undesirable local maxima.Such initial preparation is not so restrictive given that single qubit rotations are easily implemented in many quantum systems when compared with multi-qubit operations.We use values of τ = 1ns for the duration of each time interval and a total number of intervals afforded to the optimizing agent of T = 10.The results are shown in Fig. 2 where a highest fidelity of F (ρ Q , σ) = 0.9746 is achieved in ≈ 8ns.b.Three-qubit Dicke state.Dicke states embody another class of genuinely tripartite entangled states, inequivalent to GHZ states [50], which have been experimentally realised, projected onto lower dimensional entangled states and employed in open destination teleportation and telecloning [51,52].Specifically we consider the three-qubit two-excitation Dicke state, sometimes referred to as a flipped W state, defined as We assume initial state preparation of |000 Q |0 cav , and allow more fine control this time with τ = 0.5ns and T = 20.The maximum fidelity achieved in this instance is F (ρ Q , σ) max = 0.9896 within ≈ 7ns, as shown in figure 3. c.Box Cluster State.Cluster states are a class of graph states useful for measurement based quantum computing [54][55][56][57].Cluster states are represented by graphs of interconnected vertices, where if one starts with each qubit in the ground state, they are defined procedurally by applying a Hadamard gate to each qubit (represented by the vertices of the graph) then applying conditional-phase gates between qubits whose vertices are connected by edges.Here we consider a so-called Box cluster state which is defined by four vertices connected by four edges to form a square.This can be written as, where CZ i, j is the controlled-Z gate between qubits i and j and H i is the Hadamard gate on qubit i. Explicitly, the state reads: As before, assuming τ = 0.5ns and T = 20 and starting from the initial vacuum state a maximum fidelity of F (ρ Q , σ) = 0.9642 was achieved (cf.Fig. 4), this time requiring the full 10ns to achieve and maintain maximal fidelity.

B. Entanglement Detection
A central question in any state-engineering scheme is the characterization of the features of the state that has been synthesized.Within the context of our investigation, the core aspect to address is the quantification of multipartite entanglement.
The task of accurately determining the amount of entanglement in a given quantum state, which is challenging in general, is made even more difficult in multipartite settings due to the hierarchical structure of entanglement in many-body systems [4,[59][60][61][62] and the need to construct convex-roof extensions of any pure-state quantifier, when dealing with mixed states [63][64][65].
A significant tool in these endeavours is embodied by entanglement witnesses, which often offer experimentally vi-able ways of detecting (or even quantifying [66]) entanglement [58,67] that are suitable for mixed states and have been used to detect genuine multipartite entanglement (GME) already for GHZ class, W-class and graph states [53,68,69].Whilst high fidelity with a maximally entangled state is a good indicator that we have generated entanglement close to the right structure, it is useful to quantitatively check for GME in the system.A natural approach is to use so-called "fidelity based" entanglement witnesses.These are of the general form where |ψ is the state of interest and c n is the maximal overlap between |ψ and all bi-separable states.Thus, any state for which Tr(ρ ŴF ) ≥ 0 is bi-separable and consequently Tr(ρ ŴF ) < 0 indicates genuine multipartite entanglement.The overlap values c n have already been calculated for 3-qubit GHZ and W states, and 4-qubit linear cluster states [53,58], which up to local unitaries and swaps coincide exactly with the three cases considered above.In light of the analysis reported in Figs. 2, 3, and 4 these witnesses are readily implemented.Considering Tr(ρ ŴF ) = c n Tr(ρ) − Tr(ρ|ψ ψ|) < 0, where of course Tr(ρ) = 1 and Tr(ρ|ψ ψ|) = ψ|ρ|ψ is the fidelity between ρ and |ψ .This is clearly fulfilled when F (ρ, |ψ ) > c n so we can simply highlight the threshold value (c n ) of fidelity above which GME can be detected.The region where this happens is highlighted in green in each of the figures.
One may also be interested in not only ensuring that the state has GME but also, in the GHZ case, if we have generated GHZ-class entanglement [70].In this case, c n will be the maximal overlap between |GHZ and all W-class entangled states.This region is highlighted in blue in Fig. 2.
Finally, other constructions of witnesses are useful for specific states; for example, witnesses based on collective spin operators have been used to more efficiently detect GME in symmetric Dicke states as in [53].Here the witness takes the form where Ĵk is the collective k = x, y spin operator [68].It is not so straightforward here to place a delineation at a given fidelity as with the fidelity based witnesses, so we first plot Tr(ρ Ŵ) as an inset within figure 3 and highlight the region for which Tr(ρ Ŵ) < 0. This region is then shown as the grey hatched area on the larger fidelity plot.It can be seen that both the fidelity based witness and the collective spin based witness detect GME at strikingly similar times.

C. Decoherence
So far, we have exclusively considered closed system dynamics and as such the optimality of the control schemes presented is limited to the noiseless case.It is of interest then to assess how these control schemes perform in the presence of decoherence.Specifically we can write the following Lindblad master equation to model the effect of decay and dephasing acting on each of the constituent subsystems [71,72] where κ is the cavity decay rate, γ φ, j and γ j are the dephasing and decay rate for the j th qubit, respectively, and we have introduced the superoperators for an arbitrary operator Q.If we adopt physically reasonable values for these rates we can obtain an estimate of the performance of a physical system.For example, considering superconducting systems we set κ = 2π×5 kHz for the cavity damping and the typical values of 2π×300 KHz and 2π×5 MHz for the dephasing and decay rate for each of the qubits [17,23].In these systems, the use of high-Q cavities and low temperatures leads to a reduced cavity decay and dephasing rate.The qubit decay is thus the main source of decoherence.In these conditions, the effects of noise is reported in Fig. 5.We can clearly see that the control protocols are almost completely insensitive to cavity decay and qubit dephasing, whilst still reasonably robust against qubit decay.

VI. CONCLUSIONS
We have investigated the use of evolutionary algorithms, which are well-established strategies for both classical and quantum control, for direct quantum state engineering and multipartite entanglement generation.Specifically, we considered the effective Hamiltonian of a network of noninteracting qubits jointly addressed by a common driven bus and applied a genetic algorithm to identify the set of optimal pulses to drive the evolution of the qubit register.This has allowed us to successfully put forward robust protocols for the engineering of three-and four-qubit states that play a crucial role in quantum metrology and computing, including Dicke and cluster states.The protocols, which offer significant robustness to the most common and crucial sources of imperfection, provide further evidence of the benefit of a hybrid approach to quantum control that puts together the insight provided by machine learning strategies to well-established schemes for optimal control.The extension of these approaches to larger registers and non-unitary dynamics will pave the way to quantum process engineering enhanced by machine learning and optimised by quantum control methods.ing the tanh function along the y-axis and clipping it at ±1.So the time dependence within time interval i, t i ≤ t ≤ t i+1 can be written as, A1) where u i , u i+1 is the value of the control at times t i , t i+1 respectively, W determines the severity of the step, and S is a scaling factor introduced to deal with the error ε, as in fig A. Therefore the complete functional for a general control, under this scheme, is given by where, again, T is the number of time intervals.Full functional forms can be clearly seen in figures 2, 3, and 4.
Appendix B: Algorithm Implementation Details.
In each case the population size for the algorithm was N pop = 48 and was determined by the number of CPUs available, since each chromosome in the population was evaluated in parallel to significantly speed up computation time.The number of survivors was fixed at N pop /2 = 24, from which N pop /4 = 12 pairs of parents were selected according to a probability distribution determined by their relative fitness.Each pair of parent chromosomes produced 2 offspring chromosomes to re-populate.The mating procedure involved 2 steps.After making a one copy of each parent chromosomes: 1) Each separate section of T + 1 elements -corresponding to the different control pulses -were completely swapped between the two copied chromosomes with ≈ 50% probability, otherwise they were left unchanged.2) random indices were then selected and the combination given by equation 9 was applied, where β was randomly sampled from [0,1] for each separate index, in an inverse manner to the two copies.This results in 2 new offspring chromosomes that are completely complementary to one and other.Mutation was then applied by selecting chromosomes (apart from the fittest) at random and random indices within these chromosomes to replace with completely random values.The rate α determined the total number of parameter values within the entire generation that were flipped and generally assumed a value of α ≈ 0.2.

Appendix C: Simulation Details
Simulations were carried out using the Numerical Schrödinger Equation solver within the QuTiP package in python [73].Thus, for the simulations we had to approximate the harmonic mode by a d-dimensional harmonic system, however care must be taken.Since the cavity is assumed to be resonantly driven, if we use a low value of d and/or too strong a drive, then the d-dimensional harmonic system will quickly become saturated, and in fact early optimizations used this fact to produce extremely good results with low dimensional cavities.Of course when one then simulated these controls with higher levels in the cavity the performance was destroyed and so the results were neither realistic nor practically useful.One can avoid this in two ways: 1) By encoding enough redundancy in the cavity by using very high values of d, which increases computational cost; or 2) Including a term in the numerical fitness function that punishes population of the higher energy states of the cavity approximation.Here we apply both by choosing d ∈ [5,6] as well as including the term where |n is the highest excited state within the d-dimensional approximation and ρ cav is the reduced state of the cavity.ν determines the strength of punishment and was set at ν = 0.1.(Practically, since the dynamics is solved numerically the integral was actually a summation).Another issue one encounters with this type of simulation is that, if we assess fitness based on the outright maximum value of the fidelity during the induced dynamics, then it is possible to observe successful control schemes that induce sharp spikes in the fidelity landscape.If the control scheme induces such spikes on a time scales shorter than that required to completely "switch off" all of the controls then it is impossible to extract the state of maximum fidelity and the control sequences again cease to be practically useful.We can combat this by including an additional term in the fitness function that rewards control schemes that briefly maintain near maximal fidelity for a short time, allowing us to selectively uncouple the cavity whilst maintaining the state of maximal fidelity in the qubit subspace.The term used was where m is the number of time intervals of length τ to include in the numerical bonus beyond which we no longer care if the fidelity deteriorates.µ is again a variable that determines the relative importance of maintaining fidelity after maximum and was set to µ = 0.5.Thus the actual reward function employed was Chromosome fitness = F (ρ Q (t max ), σ) + φ 1 + φ 2 (C3) where t max is the time at which maximum fidelity is achieved during the induced dynamics.This leads to control schemes that maximise and briefly maintain fidelity allowing us to selective switch of the couplings, whilst also only exclusively utilising lower lying levels of the harmonic mode.Thus in principle the resulting controls could be practically implemented and yield identical performance.Clearly, these considerations are necessitated by the use of simulation and the case is much simpler if one wishes to use the algorithm on a physical system, however in this case the ability to parallelise the computational steps, one major advantage of the Genetic Algorithm, is suppressed.

FIG. 1 .
FIG. 1. Schematic representation of how the variables within each chromosome (visually represented by the sequence of connected squares) are assigned to which control pulse, depicted along the top of the diagram.Along the bottom shows sample control pulses generated using random parameter values and the function construction method outlined in Appendix A.

FIG. 2 .
FIG.2.Results for three-qubit GHZ state.We have the optimal control pulses [Top], the fidelity in the qubit subspace during the dynamics [Middle] and the matrix histogram for the target state (left) and the state with the highest fidelity during the dynamics (right)[Bottom].Hyper-parameters τ = 1ns and T = 10 are used here, with maximal fidelity achieved, and maintained, within ≈ 8ns.On the fidelity plot, the horizontal green line is drawn at c n = 1/2 and all fidelities above this exhibit GME (green region).The blue region c n = 3/4 on the other hand shows those fidelities that exhibit both GME and GHZ-class entanglement in particular.

FIG. 3 .
FIG. 3. Results for three-qubit Dicke state.We have the optimal control pulses [Top], the fidelity in the qubit subspace during the dynamics [Middle] and the matrix histogram for the target state (left) and the state with the highest fidelity during the dynamics (right) [Bottom].The hyper-parameters τ = 0.5ns and T = 20 have been used for these calculations.The corresponding maximal fidelity was achieved within ≈ 7ns.The green region on the fidelity plot highlights fidelities for which genuine multipartite entanglement (GME) is detected based on the fidelity based witness Eq. (15) with c n = 2/3.The inset shows Tr(ρ Ŵ) at each instant of time, and the grey region on both the inset plot and the fidelity plot show the temporal region where GME is detected via the collective-spin witness in Eq. (16) with b s = 3.12[53].

FIG. 4 .
FIG. 4. Results for a four-qubit Box Cluster state.We have plotted the optimal control pulses [Top panel], the fidelity in the qubit subspace during the dynamics [Middle panel] and the matrix histogram for both the target state (left-most figure) and the state with the highest fidelity during the dynamics (right-most one) [Bottom panel].We have used the hyper-parameters τ = 0.5ns and T = 20 and the maximum fidelity was achieved within ≈ 10ns.The values of fidelity for which GME is detected via the fidelity based witness Eq. (15) with c n = 1/2 [58] are shown in the green region.

FIG. 5 .
FIG. 5. A comparison between the fidelity achieved for the (Top) GHZ state preparation, (Middle) Dicke state preparation and (Bottom) Box Cluster state preparation schemes in the presence of cavity decay and qubit dephasing only (Grey) and qubit decay only (Red).Each plot shows the maximum achieved fidelity in the latter case.The horizontal blue line shows the maximum fidelity achieved in the ideal, noiseless case.

FIG. 6 .
FIG. 6.Comparison of how the clipping window affects the shape of the interconnecting tanh functions used to construct the control functions.In (a) the tanh function is shown with over-laid dashed lines representing different widths of this clipping window.From (a) we can see how the clipping window width determines the severity of the time dependence.Namely, taking the largest width (grey) leads to a more "step-like" dependence, as evidenced in (b) where each of the clipped functions are scaled into a time interval of equal duration τ and re-scaled to account for the error ε.(a) shows how the error ε increases as the clipping window narrows, i.e the narrowest (yellow) window has the largest error thus requiring the most re-scaling.