Derivative quantum circuits

Abstract

Method and systems are disclosed for computing a derivative of a quantum circuit using a hybrid quantum computing system. The method comprises receiving or determining a formulation of a parametrized quantum circuit parametrized by a continuous variable. The quantum circuit may encode a given parameter value of the continuous variable in a Hilbert space associated with a quantum elements of the quantum computer. The method further comprises determining or receiving a number of one or more distinct gap values and an equal number of distinct parameter shifts. The number of gap values is lower than the number of spectral gaps of a generator of the quantum circuit. The method further comprises executing, for each of the distinct parameter shifts, the parametrized quantum circuit with the continuous variable x increased by the respective parameter shift and with the continuous variable x decreased by the respective parameter shift. The execution of the parametrized quantum circuit comprises translating the parametrized quantum circuit into control signals for controlling the plurality of quantum elements and for readout of the plurality of quantum elements to obtain hardware measurement data; controlling the quantum computer system based on the control signals; and receiving, in response to the execution of the quantum circuit, the hardware measurement data. The method further comprises processing the hardware measurement data to obtain the derivative of the quantum circuit with respect to the given parameter.

Description

[0001]Derivative quantum circuits Technical field This disclosure relates to derivative quantum circuits, and in particular, though not exclusively, to methods and systems for determining and executing derivative quantum circuits, and to a computer program product enabling a hybrid classical quantum computer system to perform such methods. Background Quantum computing is a burgeoning field with promises of implementing algorithmsfor certain problems with a speed-up over their classical counterpart — known as quantumadvantage. These speed-ups come from properties unique to quantum computers — such assuperposition, entanglement, and an exponentially scaling computational space — whichlead to a fundamentally different way of computation. Fields which are promising for quantumadvantage to be found include fields in which quantum-mechanical phenomena play a roll, such as certain chemistry applications and molecular simulations. Other examples are provided by fields with (computationally) large problems, such as generative modelling and many differential equations. Differential equations form the basis of modelling a great number of systems, such as electromagnetic waves, evolution of heat distributions, population dynamics, optimisation of industrial processes, and time-evolution of probability distributions. Consequently, there has been much recent interest on biasing neural networks towards differential equation solutions. A common means of representing differential equations in neural networks is using Physics- Informed Neural Networks (PINNs), which use (potentially higher-order) derivatives of neural networks outputs with respect to their inputs to construct appropriate regularisation terms towards a prescribed differential equation. Calculation of the derivative of the output of a quantum circuit with respect to variable circuit parameters (input or variational parameters) is essential to implement a large class of quantum algorithms, notably quantum machine learning algorithms like Derivative Quantum Circuits (DQC) and Quantum Extremal Learning (QEL). Analytical differentiation methods like automatic differentiation are frequently favoured over numerical differentiation methods like finite difference method even in classical machine learning due to their robustness against errors. To solve differential equations or to perform gradient-descent-based optimization, it is essential to have access to the derivatives of quantum computations. Two techniques are widely used in the art: finite differencing and analytic differentiation. Finite-difference based methods compute a derivative ^^ of a function ^ at a point ^with some variation on ^^(^) ^(^^^^) ^^(^) is the finite difference. This expressionis exact in the limit that the Finite differencing is typically used ingrid-based implementations, where Δ^ is the grid step. Hence, numerical differentiationtechniques require small difference Δ^ to reduce the so-called bias error or discretisationerror, leading to a prohibitively large computational cost. However, using current or near-term hardware, it is difficult to accurately implement small parameter shifts, i.e., to accurately input Δ^. Finite difference methods may also be used in grid-based methods, where the girdprovides a natural step size. In such implementations, the computation of ^^^^ − ^^, where ^denotes a grid index, requires so-called shifter circuits to compare the function values in adjacent grid cells. However, current hardware implementations of these shifter circuits areextremely noisy; that is, the resulting error in ^^^^ − ^^ is much larger than the uncertainty in^^. Moreover, quantum measurements are inherently probabilistic, leading to anuncertainty in the determined function values ^(^ + Δ^) and ^(^). For small differences,these uncertainties are divided by a small number (the step size), and thus blow up.Consequently, there is an intrinsic limit on the accuracy of derivatives using finite difference methods for quantum computers. These uncertainties are often neglected in simulations of quantum circuits, but are a severe problem in actual hardware implementations. Thus, in general, a value for Δ^ that is too large will result in a large bias error, and avalue of Δ^ that is too small will result in a large variational error and in errors in Δ^.By contrast, analytic differentiation is, in theory, exact in the sense that the bias erroris 0. The most common way to implement analytic differentiation is by use of the so-calledParameter Shift Rule (PSR), or the Generalised Parameter Shift Rule (GPSR) as describedin Kyriienko et al., ‘Generalized quantum circuit differentiation rules’, Phys. Rev. A 104(052417). These use typically fairly large parameter shifts, so that the variational error isrelatively small (but still non-zero), and so that (hardware) errors in the exact parameter values have a limited effect. However, the (normal) parameter shift rule is only applicable to certain quantum circuits, implemented on non-interacting hardware where the parameters control single-qubit gates. The parameter shift rule cannot be used for interacting systems, such as a neutralatom quantum computer and many other current hardware implementations. The (generalised) parameter shift rule can, in principle, handle any quantum circuit, but may result in a very large number of circuit repetitions; typically at least two for everyencoding qubit, more for more expressive feature maps. In theory, for a system with ^qubits, this can require 2^^circuit evaluations (not even counting the number of shots foreach circuit to reduce the variance). This limits their practical application on currenthardware. The output in quantum machine learning, even in case of completely ideal hardware, would contain the fundamental shot noise. This, combined with the limited shot budget of near-term quantum hardware, make numerical differentiation error prone and unfeasible.. The generalized parameter shift rule provides a method to calculate the exact derivative with respect to parameters represented by an arbitrary generator with a rich spectrum of Eigenvalues while remaining robust against noise. However, the generator of aquantum circuit with ^ qubits will have 2^ eigenvalues, which corresponds to makingmeasurements at 2^(2^ − 1) different values of the parameter to estimate the derivative atany point. This quickly becomes prohibitively expensive. For example, even for a very small, five qubit circuit on neutral atom hardware, this corresponds to 992 measurements. Given a typical shot budget of about 1000 shots for a single simulation, this results in only a single measurement per measurement point, which is insufficient to determine the expectation value. There is therefore a need in the art for a method for computing (includingapproximating) a derivative of a quantum circuit with a limited shot budget, in a system with anon-trivial spectrum of eigenvalues (more than two distinct Eigenvalues). Summary It is an aim of embodiments in this disclosure to provide a system and method for computing a derivative of a quantum circuit that avoids, or at least reduces the drawbacks of the prior art. In an aspect, embodiments in this disclosure relate to a method for computing a derivative of a quantum circuit using a hybrid data processing system comprising a classical computer system and a quantum computer system, the quantum computer system comprising a plurality of quantum elements. The method comprises obtaining, by the classical computer system, a formulation of a parametrized quantum circuit parametrized bya continuous variable in a variable space. The quantum circuit may encode a givenparameter value of the continuous variable in a Hilbert space associated with the plurality of quantum elements. The method further comprises obtaining, by the classical computer system, a number of one or more distinct gap values and a number of one or more distinctparameter shifts. The number of one or more distinct gap values may be smaller than anumber of spectral gaps obtained from the Eigenvalue spectrum. The number of one or moredistinct parameter shifts may be equal to the number of one or more distinct gap values. Themethod further comprises executing, by the classical computer system, for each of the one or more distinct parameter shifts, the parametrized quantum circuit with the continuous variable shifted based on the respective parameter shift, e.g., with the continuous variable increased by the respective parameter shift and with the continuous variable decreased by therespective parameter shift. The method further comprises processing hardwaremeasurement data received from the quantum computer to obtain the derivative of theparametrized quantum circuit with respect to the given parameter. In an embodiment, the execution of the parametrized quantum circuit comprises translating the parametrized quantum circuit into control signals for controlling the plurality of quantum elements and for readout of the plurality of quantum elements to obtain hardware measurement data; controlling the quantum computer system based on the control signals; and receiving, in response to the execution of the parametrized quantum circuit, thehardware measurement data. The hardware measurement data may be representative of astate of the plurality of quantum elements. The classical computer obtaining information —e.g., the formulation of a parametrizedquantum circuit, the one or more distinct gap values, or the one or more distinct parametershifts— may comprise the classical computer receiving said information from an externalsystem (e.g., via a communication interface), or the classical computer determining said information, e.g., using software stored on a computer-readable storage medium of the classical computer. As used in this application, the derivative may refer to an exact or to an approximated derivative. As used herein, a spectral gap is defined as a positive difference between twoEigenvalues of a generator defining the dependence of the quantum circuit on the parameter.Thus, a generator with more than one distinct spectral gap may be alternatively defined as agenerator with more than two distinct (possibly degenerate) Eigenvalues. As used herein, agenerator is a Hermitian matrix that is related to a unitary operator describing a quantumcircuit via an exponential mapping. In practice, the generator corresponds to the Hamiltonianof the system, at least up to a constant multiplicative factor (conventionally, the generator is defined as two times the Hamiltonian). As used herein, a gap value is a mathematical parameter that is used in a way that is analogous to the use of spectral gaps in the (generalised) parameter shift rule. A gap value may be equal to a spectral gap of the generator of the quantum circuit under consideration, but this need not be the case. This method allows to calculate the derivative for a parametrized quantum circuitwhich can be represented by a unitary operator ^^(^) = exp^−i ^ ^^ / 2^, where ^^ is aHermitian matrix referred to as the generator of the unitary operator. Generally, the gap values are at least informed by an Eigenvalue spectrum of the generator ^^. The method strikes a balance between the bias error, variance, and number of circuit evaluations. This isachieved by using one or more parameter shifts that are based on the spectrum of thegenerator associated with the circuit. It has been found that, surprisingly, the relative error in the computed derivative increases only very slowly, or for certain ranges of number ofdistinct gap values not at all, with a decreasing number of distinct gap values. Hence,typically, the number of distinct gap values is much smaller than then number of distinct spectral gaps. In an embodiment, the method further comprises receiving or determining, by theclassical computer system, spectral information describing the Eigenvalue spectrum. In suchan embodiment, the one or more distinct gap values may be based on the spectralinformation. The spectral information may comprise, e.g., an estimate for the largest spectral gap in the Eigenvalue spectrum, a statistical description of the Eigenvalue spectrum or of an approximation thereof, a description of a dominant structure of the Eigenvalue spectrum, a complete list or substantially complete list of all spectral gaps, et cetera. Additionally or alternatively, the spectral information may comprise information about the Eigenvalues from which properties of the Eigenvalue spectrum may be derivable, e.g., estimates (or exact values) of the smallest and largest Eigenvalues, a statistical distribution of the Eigenvalues, or a complete list of all Eigenvalues, et cetera. The gap values may or may not coincide with one or more of the plurality of spectral gaps. In an embodiment, the spectral information comprises at least one of: an exact and / orapproximate value of largest spectral gap and, optionally, a smallest positive spectral gap,and / or a smallest and a largest Eigenvalue; a statistical representation of a distribution ofpositive spectral gaps; and a list of all positive spectral gaps.It has been shown that good results may be obtained when selecting gap values from the range of spectral gaps of the generator. In general, the smallest spectral gap may be assumed to be close to zero, so that just the largest spectral gap (i.e., the difference between the largest and smallest Eigenvalues) is sufficient to estimate the range of spectral gaps. More detailed spectral information allows to further optimise the selection of gap values. In some cases, smaller errors may be obtained by selecting the gap values such that a distribution of the gap values follows the distribution of the spectral gaps.In an embodiment, the one or more distinct gap values are obtainable by at least oneof: uniform or quasi-uniform sampling in a predetermined range and / or a range defined bythe spectral information; and sampling based on the spectral distribution. For example, thegap values may be sampled from a predetermined range, which may depend on the number of gap values to be selected. Additionally or alternatively, the gap values may be selected from a range defined by the spectral information, for example, a range from zero to a maximum spectral gap or a subrange thereof, or a range comprising a predetermined fraction of the spectral gaps. This range may be referred to as a range of potential spectral gaps. In an embodiment, the generator is based on a non-interacting part corresponding to a plurality of non-interacting quantum elements and an interacting part based on an interaction of the plurality of quantum elements. The contribution of the interacting part to the spectral gaps of the generator may be relatively small compared to a contribution of the non-interacting part. In such an embodiment, the spectral information may comprise informationon the Eigenvalue spectrum of the non-interacting part, and the one or more distinct gap values are obtainable by at least one of: uniform or quasi-uniform sampling in a range of potential spectral gaps of the non-interacting part; and sampling based on the distribution ofpositive spectral gaps of the non-interacting part.In many cases, the spectrum of the generator not taking into account any interaction between the quantum elements, can more easily be determined than the spectrum of the generator taking into account the interaction between the quantum elements. The interactingpart may be considered a perturbation of the non-interacting spectrum. The interactionbetween the quantum elements can be any kind of normally used interaction, e.g., Rydberg interaction, XY interaction, and so on. For example, the generator is based on a base generator corresponding to a plurality of non-interacting qubits and a contribution based on an interaction of the plurality of quantum elements, wherein a contribution of the interaction to the spectral gaps is relativelysmall compared to the base generator. This may result in clusters of spectral gaps, clusteredaround the spectral gaps of the base generator. Due to the interaction between the quantum elements, these clusters may be broadened and / or shifted. The gap values may be determined based on the spectral gaps of the base generator, possibly taking into account the effect of the interaction on the shape and central values of the clusters. In an embodiment, the distinct parameter shifts are obtainable by at least one of: selecting the parameter shifts such that the variance of the derivative is minimized; andselecting, preferably uniformly selecting, the parameter shifts from a predetermined range,^ ^^ ^ ^ preferably from the range ^ ^ , ^ ^, more preferably from the range ^ ^ , ^ ^. The parameter shifts may be based on the gap values, in particular when the variance of the computed derivative is minimized. In general, parameter shifts that are too large may increase the bias error, and parameter shifts that are either too small or too large may increase the variance. It has been found that typically, parameter shifts from the cited ranges result in an acceptable total error (where the total error depends on both the biaserror and the variance). Uniform sampling over the range is computationally straightforward,and ensures that the selected parameter shifts are not too close to each other. It has beenfound that uniform or pseudo-uniform sampling is particularly efficient (in terms of relative error as a function of number of distinct gap values) if the entire range of spectral gaps issampled, i.e., if the largest sampled value is close to a maximum of the distribution ofspectral gaps. The parameter shifts may be selected based on the hardware constraints of the quantum computer. In an embodiment, obtaining the derivative of the quantum circuit with respect to the given parameter comprises, for each of the one or more distinct parameter shifts, determining a difference between the two shifted function values. Based on the determineddifferences, the one or more distinct gap values and the one or more distinct parametershifts, a linear system of equations can be solved to determine a gap coefficient associatedwith each gap value. The derivative of the parametrized quantum circuit with respect to thegiven parameter can be computed based on the one or more distinct gap values and theassociated gap coefficients. The gap coefficients ^^ may be computed, for example, by solving a linear system ofequations of the form ^ ^wherein ^^ denotes the the ^thparameter shift ^^,and ^ is the number of distinct gap values Δ^ and distinct parameter shifts ^^.The derivative of the parametrized quantum circuit ^(^) with respect to the parameter^ parametrizing the parametrized quantum circuit may be computed by^wherein ^^ denotes the gap value Δ^, and ^ is the number ofdistinct gap values Δ^. In an embodiment, the number of distinct gap values is significantly smaller than thenumber of distinct spectral gaps. For example, the number of distinct gap values may be smaller than or equal to half the number of distinct spectral gaps; smaller than or equal to thenumber of distinct spectral gaps to a power of 1 / ^ with ^ a positive integer larger than 1, e.g.,^ = 2, 3, 4, …; or smaller than or equal to a base-2 logarithm of the number of distinct spectralgaps or a small multiple thereof. The number of distinct gap values may also be equal to orsmaller than, e.g., once or twice the number of quantum elements, or once or twice a square root of the number of quantum elements. In an embodiment, the plurality of quantum elements comprises at least eight, at leasttwelve, at least sixteen, or at least twenty quantum elements. In an embodiment, the number of distinct gap values is at most 128, at most 64, atmost 32, more at most 16, or at most 8. In an embodiment, the number of distinct gap values is at least 2, at least 4, at least8, or at least 16. The benefit of not using all distinct spectral gaps is greatest if the number of spectral gaps used is substantially smaller than the number of distinct spectral gaps. In an embodiment, the number of distinct gap values is based on an relative errorassociated with the number of distinct gap values. It has been found that the relative errorintroduced by not using all distinct spectral gaps is mostly independent of system size, anddepends only weakly on the number of gap values used. Hence, it is possible to balance theexpected error due to the approximation of not using all spectral gaps with other uncertainties, in particular uncertainties inherent in quantum systems. For example, for agiven shot budget, a balance may be struck between using a larger number of gap valueswith fewer shots per gap value (and hence a larger measurement uncertainty per gap value),and using fewer gap values and a larger number of shots per spectral gap.Although the particulars depend on the system, the number of distinct gap valuesmay be determined based on one or more of: shot budget, expected uncertainty associated with the number of distinct gap values, and system size. In an embodiment, executing the one or more parametrized quantum circuits comprises translating each of the one or more parametrized quantum circuits into a sequence of signals (e.g., pulses, typically electric or electromagnetic pulses such as voltage pulses, optical pulses, or microwave pulses) and using the sequence of signals to operate the quantum elements in the quantum register of the quantum processor. Thus, executing the one or more parametrized quantum circuits typically comprisesthe classical computer system controlling the quantum processor to perform operationsdefined by the one or more parametrized quantum circuits. Generally, at least one of the signals in the sequence of signals depends at least inpart on the given parameter. For example, a duration and / or an amplitude of the pulse maydepend on the given parameter. In an embodiment, receiving hardware measurement data comprises applying a read- out signal to quantum elements of the quantum processor and, in response to the read-out signal, measuring quantum hardware measurement data. In general, the embodiments described herein can be implemented as sequences of operations that can be executed on, e.g., a neutral atom quantum computer, either in digital mode of operation, in analog mode, or in digital-analog mode. In particular, a hybrid quantum computer as described herein may include one or more of: a gate-based qubit device, a digital / analog quantum device, a neutral-atom-based quantum device, an optical qubit device and / or a gaussian boson sampling device. As an example, the parametrized quantum circuit may be represented by a function^(^) parameterized by a variable ^. The function ^(^) is generated by a generator ^^, i.e., the^ ^-dependence of the function is of the form exp(−i ^ ^^). An approximation for the derivative^′(^) can be obtained by analysing a subset of the Eigenvalue spectrum of ^^ and thestructure of ^^. Following the Kyriienko et al., this allows the formation of a “stencil for thederivative” with unknown coefficients. This so-called stencil of the derivative refers to a sumof trigonometric functions, e.g., sine functions, of the parameter ^ of various frequencies, thatis, ∑^ ^^ sin(^^ ^ + ^^) where ^^, ^^, and ^^ are arbitrary parameters referred to as the amplitude, frequency, and This expansion in trigonometric functions (in particular, complex exponentials) isnatural in the sense that unitary quantum operators can be written as an exponentialcontaining some operator (known as the generator). By going to its Eigenbasis, a sum ofexponentials is obtained, containing spectral gaps in exponentials. Instinctively, the fact that quantum systems deal with oscillating quantum elements and complex exponentials, may serve to explain that such an expansion is natural. The unknown coefficients ^^, ^^, and ^^ in the “stencil for the derivative” forms asystem of linear equations. The the circuit at one or more shifted values of ^ may beused to find the unknown coefficients in the “stencil for the derivative” and hence thederivative itself. These shifted values can be decided, at least in part, based on the spectralgaps and the structure of ^^ obtained above. Consequently, evaluation of the quantum circuitat one or more shifted values of ^ allows determination of the derivative.The shifted values can be selected, for instance, based on the Eigenvalue spectrumof the generator, or in a uniform manner, e.g., centred around a value. Depending on theavailable information about the spectrum and / or the problem, the shifted values may beoptimized for minimum variance and / or bias.The spectral gaps (or gap values) of the methods described herein correspond to thefrequencies ^^ of this “stencil”, and can be selected based on one or more of:– uniform sampling with an upper and lower cutoff based on the spectral gaps and / orthe desired resolution in a problem the quantum circuit is trying to model; –selection of “dominant frequencies” that have a large contribution to the derivative,which may be determined depending on the encoding of the parameter ^ (this canbe done, for instance, in case of the so-called tower feature map); –selection of representative frequencies, e.g., in the case of clustered frequencies,selection of representatives of one or more clusters. If all the spectral gaps are included, and a trigonometric series is used as the “stencilfor the derivative”, this is exactly the procedure prescribed by GPSR. However, the embodiments described herein typically only use a relatively small number of gap values, which may or may not coincide with a subset of the spectral gaps; indeed, in general, it is not needed to have access to the complete Eigenvalue spectrum. Compared to finite difference methods, which have no information about thegenerator or the parameter, the embodiments described herein can obtain a lower variance for the same number of shots and a lower bias error. In particular, the embodiments described herein allow optimising the gap values based on the generator, whereas finite difference methods do not prescribe a way to optimise the step size. In general, for finite difference methods, a small step size is required to obtain an acceptable bias error, but dueto the quantum effects, a small step size leads to a high variance.The normal parameter shift rule (PSR) is state of the art for non-interacting systems,but is not compatible with interacting systems or other with generators that are notidempotent or involuntary. The generalised parameter shift rule (GPSR) works for any generator (includinginteracting systems), but requires too many circuit evaluations for generators whoseEigenvalue spectra have more than a few spectral gaps, taking into account the currenthardware limitations. Moreover, the full generator spectrum, which is a prerequisite for thegeneralised parameter shift rule, may not always be known.There are two types of relevant errors while calculating the derivative, bias andvariance. Bias refers to the error in the estimate assuming you can perform infinite shots onthe quantum hardware. Analytic differentiation methods such as the parameter shift rule have zero bias error. For finite difference methods, the bias error depends on the size of the finite difference, and tends to zero bias error only in the limit of the finite difference going to zero. Since performing infinite shots on hardware is impossible, the statistical nature ofoutput of a quantum circuit will contribute to an error in the estimate. This error is referred toas variance. The methods described herein strike a balance between shot count, bias, and variance. They may thus provide a way to compute a derivative of a quantum circuit with sufficient accuracy and acceptable computational intensity. In a further aspect, embodiments in this disclosure relate to a hybrid data processingsystem for computing a derivative of a quantum circuit. The hybrid data processing systemmay comprise a classical computer system and a quantum computer system, the classicalcomputer system comprising at least one processor and a memory storing computer- readable instructions, the memory being communicatively coupled to the at least one processor, the quantum computer system comprising a quantum register with a plurality ofquantum elements. The hybrid data processing system is configured to: determine or receivea formulation of a parametrized quantum circuit parametrized by a continuous variable in a variable space, the parametrized quantum circuit encoding a given parameter value of the continuous variable in a Hilbert space associated with the plurality of quantum elements;determine or receive a number of one or more distinct gap values; determine or receive anumber of distinct parameter shifts, the number of distinct parameter shifts being equal to thenumber of one or more distinct gap values; executing, for each of the one or more distinctparameter shifts, the parametrized quantum circuit with the continuous variable increased by the respective parameter shift and with the continuous variable decreased by the respectiveparameter shift, the execution comprising: translating the parametrized quantum circuit intocontrol signals for controlling the plurality of quantum elements and for readout of the plurality of quantum elements to obtain hardware measurement data; controlling the quantumcomputer system based on the control signals; and receiving, in response to the execution ofthe parametrized quantum circuit, the hardware measurement data; and processing the hardware measurement data to obtain the derivative of the parametrized quantum circuit withrespect to the given parameter. In particular, the at least one processor may be configured toexecute these steps in response to executing the computer-readable instructions. The hybrid data computer system may further be configured to perform any of the method steps described above. As will be appreciated by one skilled in the art, aspects of the present invention may be embodied as a system, method or computer program product. Accordingly, aspects of the present invention may take the form of an entirely hardware embodiment, an entirely software embodiment (including firmware, resident software, micro-code, etc.) or an embodiment combining software and hardware aspects that may all generally be referred to herein as a “circuit,” “module” or “system”. Functions described in this disclosure may be implemented as an algorithm executed by a microprocessor of a computer. Furthermore, aspects of the present invention may take the form of a computer program product embodied in one or more computer readable medium(s) having computer readable program code embodied, e.g., stored, thereon. Any combination of one or more computer readable medium(s) may be utilized. The computer readable medium may be a computer readable signal medium or a computer readable storage medium. A computer readable storage medium may be, for example, but not limited to, an electronic, magnetic, optical, electromagnetic, infrared, or semiconductor system, apparatus, or device, or any suitable combination of the foregoing. More specificexamples (a non- exhaustive list) of the computer readable storage medium would includethe following: an electrical connection having one or more wires, a portable computer diskette, a hard disk, a random access memory (RAM), a read-only memory (ROM), an erasable programmable read-only memory (EPROM or Flash memory), an optical fibre, a portable compact disc read-only memory (CD-ROM), an optical storage device, a magnetic storage device, or any suitable combination of the foregoing. In the context of this document, a computer readable storage medium may be any tangible medium that can contain, or store a program for use by or in connection with an instruction execution system, apparatus, or device. A computer readable signal medium may include a propagated data signal with computer readable program code embodied therein, for example, in baseband or as part of a carrier wave. Such a propagated signal may take any of a variety of forms, including, but not limited to, electro-magnetic, optical, or any suitable combination thereof. A computer readable signal medium may be any computer readable medium that is not a computer readable storage medium and that can communicate, propagate, or transport a program for use by or in connection with an instruction execution system, apparatus, or device. Program code embodied on a computer readable medium may be transmitted using any appropriate medium, including but not limited to wireless, wireline, optical fibre, cable, RF, etc., or any suitable combination of the foregoing. Computer program code for carrying out operations for aspects of the present invention may be written in any combination of one or more programming languages, including a functional or an object oriented programming language such as Java, Scala, C++, Python or the like and conventional procedural programming languages, such as the “C” programming language or similar programming languages. The program code may execute entirely on the user's computer, partly on the user's computer, as a stand-alone software package, partly on the user's computer and partly on a remote computer, or entirely on the remote computer, server or virtualized server. In the latter scenario, the remote computer may be connected to the user's computer through any type of network, including a local area network (LAN) or a wide area network (WAN), or the connection may be made to an external computer (for example, through the Internet using an Internet Service Provider). Aspects of the present invention are described below with reference to flowchart illustrations and / or block diagrams of methods, apparatus (systems), and computer program products according to embodiments of the invention. It will be understood that each block of the flowchart illustrations and / or block diagrams, and combinations of blocks in the flowchart illustrations and / or block diagrams, can be implemented by computer program instructions. These computer program instructions may be provided to a processor, in particular a microprocessor or central processing unit (CPU), or graphics processing unit (GPU), of a general purpose computer, special purpose computer, or other programmable data processing apparatus to produce a machine, such that the instructions, which execute via the processor of the computer, other programmable data processing apparatus, or other devices create means for implementing the functions / acts specified in the flowchart and / or block diagram block or blocks. These computer program instructions may also be stored in a computer readable medium that can direct a computer, other programmable data processing apparatus, or other devices to function in a particular manner, such that the instructions stored in the computer readable medium produce an article of manufacture including instructions which implement the function / act specified in the flowchart and / or block diagram block or blocks. The computer program instructions may also be loaded onto a computer, other programmable data processing apparatus, or other devices to cause a series of operational steps to be performed on the computer, other programmable apparatus or other devices to produce a computer implemented process such that the instructions which execute on the computer or other programmable apparatus provide processes for implementing the functions / acts specified in the flowchart and / or block diagram block or blocks. The flowchart and block diagrams in the figures illustrate the architecture, functionality, and operation of possible implementations of systems, methods and computer program products according to various embodiments of the present invention. In this regard, each block in the flowchart or block diagrams may represent a module, segment, or portion of code, which comprises one or more executable instructions for implementing the specified logical function(s). It should also be noted that, in some alternative implementations, the functions noted in the blocks may occur out of the order noted in the figures. For example, two blocks shown in succession may, in fact, be executed substantially concurrently, or the blocks may sometimes be executed in the reverse order, depending upon the functionalityinvolved. It will also be noted that each block of the block diagrams and / or flowchartillustrations, and combinations of blocks in the block diagrams and / or flowchart illustrations, can be implemented by special purpose hardware-based systems that perform the specified functions or acts, or combinations of special purpose hardware and computer instructions. Brief of the drawings The embodiments will be further illustrated with reference to the attached schematic drawings, in which: Fig. 1A schematically depicts a quantum machine learning model, Fig. 1B shows anexample of a basic form of a quantum feature map, and Fig.1C shows an example of aderivative thereof according to an embodiment;Fig. 2 is a flowchart of method according to an embodiment;Fig. 3 is a flowchart of part of a method according to an embodiment;Fig. 4 is a flowchart of another part of a method according to an embodiment;Fig. 5A and 5B are flowcharts of yet other parts of a method according to anembodiment; Fig. 6A and 6B depict a computer system for executing methods for computing aderivative according to an embodiment of the invention; Fig. 7 is a graph showing an expected accuracy of a method according to anembodiment; Fig. 8A and 8B are graphs showing exemplary spectra;Fig. 9A–9D are graphs showing experimental results of various methods according toembodiments; Fig. 10A and 10B are graphs showing exemplary spectra;Fig. 11A and 11B are graphs showing experimental results of various methodsaccording to embodiments; Fig. 12A–12D are graphs comparing a method according to an embodiment with aprior art method; Fig. 13 is a graph comparing a method according to an embodiment with a prior artmethod; Fig. 14A is a graph showing an exemplary spectrum and Fig.14B is a graph showingexperimental results of various methods according to embodiments;Fig. 15A–E schematically depict hardware-level schematics illustrating the applicationof logical operations to qubits using a quantum circuit; Fig. 16A and 16B schematically depicts a quantum circuit based on a digital-analogimplementation; Fig. 16C schematically illustrates pulses corresponding to operation of aquantum circuit; Fig. 17A and 17B schematically depicts a quantum circuit based on a further digital-analog implementation; Fig. 18A and 18B schematically depicts a quantum circuit based on another digital-analog implementation;Fig. 19 schematically depicts a quantum circuit based on a yet another digital-analogimplementation; Fig. 20A and 20B are block diagrams illustrating, respectively, an exemplary hybriddata processing system and an exemplary classical data processing system that may be used for executing methods and software products described in this disclosure; Fig. 21A and 21B are graphs showing a required number of gaps and a relative error,respectively, as a function of system size; andFig. 22A and 22B are graphs showing a required number of gaps and a relative error,respectively, as a function of system size. Detailed description Fig. 1A schematically depicts a known method to solve a differential equation. Theshown method is based on variational solving the problem using differentiable quantum circuits (DQC). As shown in the figure, the quantum model for the parametrised targetfunction ^^(^) may be represented by a quantum circuit including a feature map ^^^(^) 1041for encoding an independent variable ^ into the Hilbert space of a quantum register 102, andvariational Ansätze ^^^^(^) and ^^^^(^) 106, parameterized by trainable parameters ^. Thus,one can write: =^^ ^ ^^ ^ ^^ ^ ^^where ^^ The derivative of the target function is obtained using the same variational Ansätzewith the same values for trainable parameters ^ but with shifted values of ^ (the exactmethod to obtain this derivative is described later). A solution to the differential equation is found by variationally optimizing the trainable parameters ^. The quantum circuits may include sequences (layers) of quantum gates which may be executed on the quantum register. Here, the only dependence on the input value (independent variable) is in the feature map; the variational Ansatz is independent of the input value ^; and conversely, the feature map is independent of the trainable parameters ^. The training of the quantum model may be based on an optimization scheme, wherein after initialisation of the variational parameters, the gate operations of the quantum circuit may be executed, and the state of the quantum register may be measured using acost function. Based on measured states 1081,2, a loss value may be computed based on aloss function (block 110). The loss function is based on the differential equation. For example, if the differential equation reads: ^^^(^) − ^ ^(^) − ^(^) = 0, (2)(wherein ^(^) is a known function of ^), then the loss function contains a contribution basedon measurements 1081 representing ^(^) and a contribution based on measurements 1082representing ^^(^). The scaling factors ^ and ^ and the known function ^(^) can be includedby the (classical) optimiser without any quantum measurements. As the feature map encodes only a single ^-value at a time, the loss function comprises a sum over a number of collocation points ^^. For the above example differential equation, the loss function may be: ^^^^ ^ ^where ‖⋅‖ denotes a contribution for a boundary condition, which is omitted here for simplicity.) If the loss value does not meet one or more optimization conditions, the variational parameters may be updated by an optimizer 112. Thereafter, the quantum circuit may be executed using the updated variational parameters and the state of the quantum registered may be measured and a new loss value may be computed. These steps may be repeated until the one or more optimization conditions are met. The state of the optimized quantum circuit may represent an optimized approximation of a desired target function. Thus, the logic gates of the quantum circuits used in the quantum machine learningscheme are parameterized by a set of trainable parameters ^ which are trained by a classicaloptimizer to minimize a loss function. In general, these parametrized quantum circuits can be used as quantum machine learning models, applied to applications such as learning or generating an unknown distribution from a set of input features or solving differential equations. Fig. 1B shows an example of a basic form of a quantum feature map, which is hereillustrated as an example of a ‘time encoding’ type feature map. In this example, the featuremap 1041 is implemented using a Hamiltonian ℋ^ acting on all the qubits simultaneously andis parametrized by a function of variable ^ which translates to the time for which theHamiltonian ℋ^ is applied. Throughout this applications, qubits are used as standard exampleof quantum elements, but most arguments are equally applicable to systems with otherquantum elements, e.g., qudits. Measurement of the quantum state requires onemeasurement per encoding qubit. Depending on the quantum hardware, the qubits may be measured one-by-one, or all in a single measurement. Fig. 1C illustrates an example of a derivative quantum circuit for the time encodingfeature map of Fig.1B. Differentiation over variable ^ involves calculating the output of thecircuit with shifted values.Fig. 2 is a flowchart of a method for computing a derivative of a quantum circuitaccording to an embodiment. The method may be performed by a hybrid data processingsystem comprising a classical computer system and a quantum computer system comprising a plurality of quantum elements. For example a hybrid data processing system as described below with reference to Fig.20. In the example shown in Fig.2, the problem to be solved is to find an (exact orapproximate) derivative ^^ ^(^) of a function ^(^) represented by a parametrized quantumcircuit ^= ^^^^^^ ^^ ^^ ^^^where ^^ ^) ^^ (^ a a^^ with more Eigenvalue gaps. ^^ denotes an (^-independent) cost operatorand |^^ denotes a quantum state (also ^-independent). It is noted that for an arbitrary ^^ andeq. (4) represents a completely general quantum circuit, as any (^-independent)before or after ^^(^) may be absorbed into |^^ and ^^, respectively. As noted, in thiscontext, ^ represents the variable with respect to which the derivative is to be determined. Inthe context of a quantum neural network as described above with reference to Fig.1A, thiscan be either an input parameter ^ or one of the model parameters ^.Thus, in a first step 202, the method comprises receiving or determining a formulationof the parametrized quantum circuit to be differentiated, and an indication of the parameterwith respect to which the parametrized quantum circuit is to be differentiated, whichparameter is here represented by a continuous variable ^ in a variable space ^. Theparametrized quantum circuit may encode a given parameter value of the continuousvariable ^ in a Hilbert space associated with the plurality of quantum elements.A step 204 comprises determining or receiving spectral information describing anEigenvalue spectrum of a generator ^^ associated with the parametrized quantum circuit. Ifthe Eigenvalue spectrum comprises only a single spectral gap (i.e., of the generator has only two distinct Eigenvalues, typically ±^), the normal Parameter Shift Rule may be applied.Hence, it is assumed that the Eigenvalue spectrum comprises a plurality of distinct spectralgaps. A spectral gap is defined as a positive difference between two distinct Eigenvalues ofthe generator. For example, if the generator has three (possibly degenerate) Eigenvalues^^ < ^^ < ^^, then the potential spectral gaps are ^^ − ^^, ^^ − ^^, and ^^ − ^^ (the first twodifferences may be the same gap, e.g., for an Eigenvalue set {−^, 0, ^}). ageneratorwith more than one distinct spectral gap is equivalent to a generator with more than twodistinct (possibly degenerate) Eigenvalues. Some examples of Eigenvalue spectra, and waysto describe those spectra with spectral information, are discussed in more detail below, withreference to Fig. 10A and 10B.The spectral information may comprise, e.g., an estimate for the largest spectral gap in the Eigenvalue spectrum, a statistical description of the Eigenvalue spectrum or of an approximation thereof, a description of a dominant structure of the Eigenvalue spectrum, a complete list or substantially complete list of all spectral gaps, et cetera. Additionally or alternatively, the spectral information may comprise information about the Eigenvalues from which properties of the Eigenvalue spectrum may be derivable, e.g., estimates (or exact values) of the smallest and largest Eigenvalues, a statistical distribution of the Eigenvalues, or a complete list of all Eigenvalues, et cetera. In some embodiments, step 204 may be omitted, and the gap values may bedetermined in other ways, e.g., heuristically. Astep 208 comprises determining or receiving a number ^ of one or more distinct^ ,, (typically) based on the spectral informat ^^ gap values {Δ^}^^^ ion. The gap values ^Δ^^^^^may or may not coincide with one or more of the plurality of spectral gaps{Δ^}^^^^of the generator ^^. The manner of selecting the gap values can depend on the spectral information. For example, if only a spectral gap range is known (e.g., because the spectral gapinformation comprises only a maximum spectral gap), the gap values may be selecteduniformly or quasi-uniformly from that range. If a full range of spectral gaps is known, the selected values may be optimized numerically, as is briefly described below with reference to Fig.4. Astep 206 comprises determining or receiving one or more distinct parameter shifts{^ ^^}^^^ . These should be one parameter shift for each gap value, so typically, the number of(distinct) parameter shifts is equal to the number of distinct spectral gaps. As an example, ^ ^^ ^ ^ ^ from the range ^ ^ , ^ ^. bias error and variance. The parameter shifts may be determined based on the gap values, or independently thereof. In the following, parameter shift, spectral gaps, and gap values are assumed to bedistinct, unless explicitly or implicitly defined otherwise. Astep 210 comprises receiving information defining the quantum circuit ^(^). Thisinformation may have the form of, or be convertible to, a set of control values to implementthe quantum circuit on the quantum computer. A step 212 comprises receiving the parametershifts {^^^}^^^ . A step 214 comprises using the quantum computer to compute shifted valuesthat correspond essentially to {^(^ ± ^^)}^^^^. In some implementations, these steps may becombined into, e.g., a set of instructions for implementing the quantum circuit for ^(^ + ^^)and ^(^ − ^^).In other words, for each of the one or more distinct parameter shifts, the parametrizedquantum circuit is executed with the continuous variable ^ increased by the respectiveparameter shift and with the continuous variable ^ decreased by the respective parameter shift. Typically, for each shifted value, the measurement is repeated a number of times to obtain an expectation value with a known, or at least estimable, uncertainty. The execution of the parametrized quantum circuit comprises translating the parametrized quantum circuit into control signals for controlling the plurality of quantum elements and for readout of the plurality of quantum elements to obtain hardware measurement data, and controlling the quantum computer system based on the control signals. This results in the generation of measurement data, which can be converted to numerical values representing the output of the quantum circuit, i.e., numerical values for^(^ ± ^^).Based on the determined values {^(^ ± ^^)}^^^^, in a step 216, a set of centraldifferences {^^(^) = ^(^ + ^^) − ^(^ − ^^)}^^^^is determined. Based on these differences ^^}^^ {^^^, the set of gap values ^^Δ^^^^^, and the set of parameter shifts {^^}^^^^, in a step 218, a system of linear equations may be defined and solved for a set of gap coefficients {^^}^^^^. A typical example of such a linear system of equations is: ^ ^ Using the determined set of gap coefficients{^^} ^^^^and the set of gap values ^^Δ^^^^^, the derivative may be computed as: ^ ^ Following Kyriienko et al., a expectation of a Hermitian cost operator ^^: ^= ^^^ ^^ ^^ ^(^) is then given by^ (8)where ^ is the number of distinct Δ^ of generator ^^. Spectral gapsequal to zero do not contribute to the derivative, and may be ignored. In the following, spectral gaps will be assumed to be positive, even if not labelled explicitly as such. The gap coefficients ^^necessary for derivative evaluation can be obtained by solving the following system of linear equations: ^ (9)Here, ^^ = ^(^ + ^^) − ^(^ − ^^) between shifted functions^(^) with parameter shifts ^^.In a more compact way, the linear equation system can be written as a matrix equation: ^^×^ = ^^×^^^×^, (10)wherein ^^×^ is an ^ × ^ matrix with elements ^ ^×^^^ = ^^^^^ ^ and ^ and ^ are vectorscontaining elements ^ ^and ^^, respectively.Surprisingly, it has been found that good results may be obtained using only ^ ≪ ^terms (i.e., using ^ distinct gap values ^ Δ^ and ^ distinct parameter shifts ^^ with ^ = 1, … , ^).Mathematically, this can be understood by redefining the gap coefficients ^^to ^^^on theright-hand side of system of equations (9) and introduce new gap values ^^Δ^ ^^^^^ . In thegeneral case, the gap values^Δ^do not (necessarily) coincide with the ‘original’ set of generator spectral gaps{Δ^^}^^^ . With this notation the ^-th equation of system (9) reads:^ ^Δ^ ^^^^ ^^^(11) 2 or in the matrix notation, The leads to an approximation of derivative, given by: ^ ^ Using matrix notation, the “effective” ^^×^ = ^^×^^^×^, (13)The general solutions for the ^^^of this equation can be obtained using Cramer’s rule: det^×^^ ^=^^(^) (14)det(^^×^) ^^×^is a matrix formed by replacing ^-th column of the ^^×^matrix with the vector . can be simplified by using the explicit components of the vector ^^×^with original set of spectral gaps {Δ^}^^^^ given by eq. (9). For example, for ^ = 2, the value for ^^^can be obtained by ^^ ^^ ^where |⋅| . comes from the multilinearity property of determinants stating that if a column of matrix can be written as a sum of column vectors then the determinant of such matrix is equal to the sum of determinants with the corresponding column substituted with each of the vectors.Note that each function ^ ^ ^ ^ Δ^,^ Δ^, ^^, ^^^ only depends on the correspondingspectral gap Δ . For ^= 2 and some 1 ≤ ^ ≤ ^ the explicit form of^ ^ ^ ^^ ^^ ^^ ^^^ == = .The results of eqs. (14) and (15) can be straightforwardly generalized to givesolutions for any system containing any number ^ of equations:^ ^ det^^^(^×^^)^ (^^×(17) det^)Herein, the ^-th ^^×^composed of elements ^^×^ ^^ ^^^^^aning that ^^×^ =4 sin ^ ^, me = ∑^ ^^^ ^^(^×)^^^. Note that the functions^ ^ ^ ^^ ^^ ^ above are given by ratios of determinants, thus their explicit sine functions similarly to eq. (15). Plugging the expression of eq.(17) into the definition given in eq. (12), the approximate derivative can be written as^ ^ ^ ^ wherein the ^ To ,^products ^^ ≡ ^ ^^^, with ^ being a small parameter. It may be noted that eq. (19) isstructurally similar to eq. (8), with the functions ^^taking the place of the spectral gaps Δ^. Since ^ is a function of parameter ^, the ^^can be expanded with respect to ^. for ^ = 2, the combination of eqs. (15) (19) yields the following expansion:^^^^^ ^ This expression^Thus, if ^ is sufficiently and, consequently, the parameter shifts ^^ and ^^ are small aswell, it is sufficient to use only ^ = 2 equations to estimate the derivative while introducingonly 4thorder errors. However, as ^ gets bigger, ^(^^) and higher order terms might get significant and athird or even more equations may be needed to estimate the derivative more accurately.Repeating the same expansion procedure for ^ = 1 and ^ = 3 results in, respectively:^ = 1: ^ = ^^ Δ^ + ^ ^Indeed, with more equations − ^smaller, more specifically Δ ^^^ − ^^ = ^(^ ). Thus, eq. (12) can be written as:^ ^ This equation contains the dependence on the set of arbitrary gap values ^^Δ^^^^^only in the higher-order error terms of approximation of spectral gaps Δ^. Disregarding theseerror terms, the expression in eq. (23) coincides with the exact expression of derivative in full^ ^ GPSR form as given in eq. (8). This implies that, in principle, any gap values ^Δ^^^^^may beused to solve the linear equation system of eq. (11) whilst still obtaining an accurateapproximation of the real derivative as given in eq. (8). Terms involving ^ appear only in higher-order terms for large values of ^ and thus,these error terms are sufficiently small for larger values of ^ and, hence, for larger parametershifts ^^as well. Consequently, the more equations are included in the linear system, the larger the parameter shift that can be used to retain a similar accuracy of the derivative estimation. This follows also from application of the mathematical induction principle. Since it hasbeen shown that ^^ ≈ Δ^ for low numbers ^ of equations, it can be assumed that suchrelation holds for ^ = ^ − 1, so that the error is ^^^^(^^^)^. Therefore, for ^ = ^ the errorshould be on the order of, or smaller than, ^(^^^). Since in this case, the full equationsystem is solved, it follows that the prefactor of each ^^ is exactly equal to Δ^, or equivalently,that the error is infinitely small, i.e., ^(^^) with ^ → ∞. Hence, the expansion is valid for ^ =^, and therefore for ^ = ^ − 1 as well. Consequently, the assumed relation of errormagnitude holds for any ^. For general ^, eq. (20) can be rewritten as:^ ^where ^ is a of the form ^^Δ^^ − Δ^^^. This means that if a gap value^Δ^is selected equal to some spectral gap Δ^of the generator, the corresponding error term vanishes. This observation naturallyleads to conclusion that when ^ = ^ equations are included in the linear system of (9), and^ the gap values are selected to coincide with the actual set of spectral gaps {Δ^}^^^^, the original as given by GPSR is recovered, since all the error terms vanish. Fig. 3 illustrates various options of determining or receiving information describing anEigenvalue spectrum. Depending on the embodiment, step 204 as described above withreference to Fig.2 may comprise determining or receiving (in a step 3041) one or more representative values of the spectral gaps defined by the Eigenvalue spectrum associatedwith the quantum circuit to be differentiated. Generally, the smallest spectral gap mayassumed to be zero or close to zero, and consequently an estimate for the range of spectral gaps may be provided by (an estimate of) the largest spectral gap. In other embodiments, instead of the maximum, one or more different representative value may be provided, such as the 99thpercentile, 95thpercentile, 90thpercentile, mean, median, interquartile range, et cetera. Additionally or alternatively, the information describing the Eigenvalue spectrum may comprise a (statistical description) of a distribution of spectral gaps (step 3042). For example, the spectral gap distribution may be described, at least approximately, by a known statistical distribution. As a further example, the generator may be separable in an unperturbed generator (typically with a relatively simple spectrum) and a perturbation. This can, for instance, be the case in a neutral-atom implementation, where the unperturbed generator represents the Hamiltonian of a non-interacting system, and the perturbation represents the deviation from the unperturbed generator caused by the interaction between the neutral atoms (qubits). It is noted that such interaction is required for the quantum system to obtain quantum supremacyover a classical system. Examples of other hardware implementations are described belowwith reference to Fig.15A–E. Consequently, the spectral distribution may be based on a non-interacting part corresponding to a plurality of non-interacting quantum elements and an interacting part based on an interaction of the plurality of quantum elements, wherein a contribution of the interacting part to the spectral gaps of the generator is relatively small compared to acontribution of non-interacting part. In such a case, the spectral information may compriseinformation on the Eigenvalue spectrum of the non-interacting part. An example is describedin more detail below with reference to Fig.10A and 10B.In such a case, the gap values may be selected based on the distribution of thespectral gaps of the non-interacting part. Alternatively, information about the non-interactingpart may be used to estimate properties of the spectral distribution, such as a range of spectral values. The gap values may then be selected, for instance, via uniform or quasi- uniform sampling in a range of potential spectral gaps of the non-interacting part. Step 3043 describes a case where the spectral information comprises the full spectrum of spectral gaps. In practice, this is often not obtainable for a reasonably sized system with a non-trivial generator, but it can be useful for, for instance so-called toy systemsor systems with a relatively easy to analyse generator.Fig. 4 illustrates various options of determining or receiving gap values. The optionsfor selecting gap values may depend on the spectral information as received or determinedin step 204. In general, the more detailed the spectral information, the more options there areavailable for selecting gap values. For instance, a method that requires minimal information on the spectral distribution isuniform a quasi-uniform sampling from a range (step 4081). Examples are selecting gapvalues with fixed intervals (e.g., ^ Δ^ = ^ + ^ ^, ^ = 1, 2, … , with ^ and ^ positive constants;exemplary choices are ^ ≪ 1 and ^ ≈ 1), random from fixed intervals, or samplingaccording to a uniform distribution. Other methods are known in the art. If a more detailed distribution is known, the gap values may be selected based on the spectral distribution (step 4082). For example, the gap values may be selected such that the distribution of the gap values is similar to the spectral distribution. That way, it can be ensured the full range of spectral gaps is represented by the selected gap values. Moreover, regions with large numbers of spectral gaps may have more gap values selected. In particular, when the spectral distribution features a number of peaks, the gap values may be selected such that there is one gap value for each peak, or for a subset of peaks. This can be the case, for instance, when the generator can be written as an unperturbed generator with a relatively small perturbation; in that case, the gap values may be selected based on the spectrum of the non-interacting system (step 4083).However, the form of error terms in eq. (24) also suggests that also when thederivative is approximated by using only ^ < ^ gap values (and hence, equations in the^ linear system (11)), the set of gap values can be selected such that the error terms ^ ^ are minimized. Intuitively, this happens when the gap values are close to as many actual spectral gaps as possible, thus reducing the terms and pushing the sum in parenthesis multiplied by corresponding ^^closer to Δ^. In other words, by selecting ^ ^ the gap values ^ such that they are a good representation of the distribution of the spectral gaps^^^, the error of the approximation can be kept small. Inspired by the above results, a ^-dependent error function may be defined as: ^^(Δ^) = ^^(Δ^) − Δ^ , (25)or ^^(Δ) = ^(Δ) − Δ, dropping the subscript ^. It turns out that for many generators, thereexists an interval ^ of generator spectral gaps where ^^(Δ^) is small. Exemplary behaviourof this function ^^(Δ) is shown in Fig.7. In such cases, a good choice for the gap values^^ is a set values sampled from a uniform distribution over the interval [^, ^ − 1 + ^], where ^ is a small lower bound. For example, the gap values ^^ may be selected as the^ − 1 smallest non-negative integers, offset by a small ^, i.e., {^, 1 + ^, 2 +^, … , ^ − 1 + ^}. The gap values may also be scaled based on the range of spectral gaps,such that at least most of the range is covered. More specific examples are discussed below with reference to Fig. 10 and 11.Fig. 5A illustrates various options of determining or receiving parameter shifts. Forexample, the parameter shifts may be selected from a predetermined range (step 5062), e.g., ^ ^^ or f^ ^^ ^rom the range ^^ , ^ ^. The parameter shifts may be selected, for quasi-uniformly from the predetermined range. Alternatively, the parameter shifts may be selected to minimize the variance (step 5061). Fig. 5B is a flowchart of a method to minimize the variance of the derivative. Inparticular, Fig.5B is a flowchart of a method to select a set of ^ parameter shifts {^^}^^^^that minimize the variance of a derivative of a quantum circuit computed according to anembodiment, based on a set of ^ gap values ^^Δ^ ^^^^^. Application of quantum circuits, including those according to embodiments, on real quantum devices introduces an additional aspect affecting the accuracy of estimated derivatives—the probabilistic nature of quantum mechanical measurement processes usedto calculate the expectation value of the cost operator ^^, i.e., the function value ^(^) asdefined in eq. (4). For a finite number of measurements ^^^^^^, each measurement results ina random value of ^(^) that can be assumed to follow a normal probability distribution^ characterized by a mean value ^(̅^) and a variance^^(^)^^^^^^. The corresponding derivative^^̅(^)^^ depends not only on the number of shots, but alsoon the number ^ of equations included in the linear system of eq. (9) and on the chosen shiftvalues{^^^}^^^ , since it depends on the values of shifted functions ^(^ ± ^^). It turns out thatit is possible to minimize the variance of the derivative estimation by selecting optimal shift values{^^} ^^^^. This can be seen by first considering the full linear equation system of eq. (9)representing exact GPSR, and rewriting the definition of the derivative^^(^) ^^ to emphasize the dependence on shifted function values: ^ In order to calculate the for ^^(^^, ^^, … , ^^) are needed. Above, the linear system of eq. (9) was solved with Cramer’srule; however, for the current derivation, it is easier to consider matrix-vector equation (10) and to write the solution for ^^×^as ^^×^ = (^^×^)^^ ^^×^. (27)This expression makes clear that ^^(^^, ^^, … , ^^) linearly depends on all thedifferences ^^, so that it can be written in the following way: ^ ^ As evaluated with some measurement protocol on a real quantum device, and therefore,^^(^^, ^^, … , ^^) in eq. (28) can be interpreted as a sum of scaled independent randomvalues. Following Kyriienko et al., it can be assumed that the stochastic values ^(^ ± ^^) areidentically distributed and that the dependence on the argument value is not strong. This implies that their variances are constant: ^^(^) ≈ ^^^. Thus, the variance of ^ , which issum of independent and identically distributed variables in eq. (28), as a sumof scaled variances of individual stochastic variables: ^ ^ The derivative is given by eq. ,(29) we can finally write the variance of the derivative ^^^in the following way: ^ ^^ ^ ^ Since the coefficientsthey are functions of all the spectral gaps {Δ^} and parameter shifts {^^}. This means that wecan minimize the function ^ ^with respect to the shifts of parameter shift values {^ ^^^ ^} that results in the smallest variance of the derivative given byeq. (26) when estimated on a real quantum device. This optimization procedure can also beapplied with methods according to embodiments, when considering only a part of equations in system of eq. (9). For example, the spectral gaps are chosen as described above withreference to Fig. 4 or as described below with reference to Fig. 10 and 11, e.g., to minimizethe error resulting from the truncated equation system. Based on the selected (optionallyoptimised) gap values, the shift minimization can be performed as described above.Thus, the procedure of selecting optimal, variance-minimizing parameter shift valuesaccording to an embodiment can be summarized in the following way, as shown in Fig.5B.An optional step 512 comprises finding a set of optimal gap values that minimize ^^Δ^^the bias error in the derivative for a given number ^ of gap values (and hence, includedequations in system (9)), for example as described above with reference to Fig.4.Alternatively, a set of gap values selected in a different way may be received. Astep 514 comprises calculating an analytical inverse of the matrix ^^×^ to obtaincoefficients ^^^^^^Δ^^, {^^}^, where ^^×^ is a ^ × ^ matrix with coefficients ^^^ =^^ 516 comprises finding the set of optimal parameter shifts {^^} by minimizing aloss function ^^^^Δ^^, {^^}^ = ∑^ ^^^ ∑^ ^^^^ Δ^^^^^^ ^^^Δ^^, {^^}^(with respect to parameter shifts {^^}). An optional step 518 comprises using the received (possibly optimal) gap values andthe optimal parameter shifts to compute the derivative on a real quantum device with minimal variance and, optionally, minimal bias. Fig. 6A and 6B schematically depict a method for computing a derivative of afunction represented by a quantum circuit according to an embodiment. In particular, Fig.6Ashows a high-level schematic of the interplay in a hybrid computer system 600 between aclassical computer 602 and a quantum computer 604, while Fig. 6B depicts a schematic ofthe quantum computer 604. The quantum computer executes gate sequences called ‘quantum circuits’ 640, for example those described below with reference to Fig.16A–19. These quantum circuits typically comprise a qubit initialization step 642, followed bysequences 644,646 of digital 648,652 and / or analog gates 650. These gates can be, forexample, unitary operations on a quantum state or measurement operations 654. Examples of hardware-level schematics to implement these quantum circuits are discussed in more detail below with reference to Fig.15A–E. In the methods described in this disclosure, these quantum circuits can be used to represent a function and a derivative of that function, represented by a parametric quantum circuit, e.g., a quantum neural network, or a quantum kernel. In the depicted example, gatesequences 606 are used to encode an input parameter ^ ∈ ^ into a Hilbert space associatedwith the quantum computer. For example, for a quantum computer with ^ spin-½ qubits, thiscan be In embodiments where ^ is an independent variable or input variable,gate may be referred to as a quantum feature map circuit. However, inprinciple there is no difference, insofar as computing a derivative is concerned, with model parameters ^; for clarity, in this disclosure, the parameter with respect to which the derivative is to be computed is generally denoted ^, whereas any other parameters that may be used otparametrize the quantum circuit are denoted by ^. Additionally, e.g., subsequently orconcurrently, gate sequences 608 are applied which are parametrised by parameters ^. Inthe case of a machine learning model, the parameters ^ may be initialised with properlychosen initial values, e.g., random angles and subsequently optimised. Several rounds ofapplying gate sequences 606,608 and measurement 610 may be required to determine asingle expectation value with sufficient accuracy. The measurement may be based on (anexpectation value of) a cost function ^^. The result is output data 626.The output data 626 are used by the classical computer 602 to compute 612 afunction value ^(^). The function value may be equal to the computed expectation value, forexample, or be based on the output data in a different way; for example, via scaling, translation, or using other mathematical operators. It is noted that in this figure, thecomputation of the function value ^(^) is only shown for the sake of comparison; in general,there is no need to determine the function value in order to determine the derivative function value. The quantum circuit 616 that is referred to as a derivative quantum circuit, is generallythe same quantum circuit as described above to compute the function value, with onedifference. Instead of parametrizing the quantum circuit with the parameter ^, the quantumcircuit is parametrized with a shifted value of the parameter to be differentiated, and executed (at least) twice, once with a shift in a positive direction and once with a shift in anegative direction. This is repeated for each of a number ^ of one or more parameter shifts {^^}^^^^. Any operations parametrized by other parameters ^, e.g. a second parametrizedquantum circuit 608, remain unchanged. Again, several measurements 620 may be used toobtain an expectation value for each shifted parameter value. The number of circuitrepetitions (shots) may the same for the function value and the derivative function value, ordifferent. The result is output data 620. Based on these output data, the classical computer may determine the value for the derivative function 622. This value may be used, e.g., as input in a loss function, to optimize model parameters (e.g., using gradient descent / ascent), et cetera. Data representing the measurements 626,630 are transferred to the classical computer 602. Typically, the data are transferred only once all calculations on the quantum computer are finished. Then, the derivative function value (and optionally the function value)is determined, based on the data received from the quantum computer 604.Fig. 7 is a graph showing an expected accuracy of a method according to anembodiment. In particular, Fig.7 shows graphs for the estimated error ^^(Δ) on the verticalaxis as function of spectral gap Δ on the horizontal axis, where ^^(Δ) is as defined in eq. (25)using equidistant gap values ^ Δ^ ≈ ^, ^ = 1, 2 … , ^. Graphs are shown for ^ = 2 (graph 702),^ = 4 (graph 704), ^ = 8 (graph 706), ^ = 12 (graph 708), and ^ = 20 (graph 710).In Fig. 7, it can be seen that for all values of linear equation system (11) of size ^there is indeed an interval of gap values Δ where the error function ^^ is small (substantiallynear zero), i.e., ^ ≪ 1, e.g., (depending o −2 −4 −8^ n ^) < 10 , < 10 , or < 10 . With larger values of^, the interval where ^^ ≪ 1 grows. Fig.7 was generated using shift parameters valued ^^ ∈^ ^ ^^ that are sufficiently large.An immediate consequence of the behaviour of ^^(Δ) shown in Fig.7 is that if somegenerator ^^ produces spectral gaps described by distribution ^(^) (mostly) confined to theinterval where ^^(Δ) ≪ 1, the estimated derivative will be practically as accurate as whenusing the full GPSR method with all spectral gaps included in the solution of the linearsystem of eq. (9). This can be seen in Figs.8 and 10.Figs.8–11 show results for calculating the derivative of function defined in eq. (7)with a typical neutral-atom generator (setting ℏ = 1) of ^ qubits, given by:^ ^^^ In the generator above, atoms, ^^^^^^^^^^ is the interaction strength between atoms ^ ^ and ^^ is the Pauli operator.Fig. 8A shows a distribution of spectral gaps ^(^) for a neutral-atom systemgenerator dominated by inter-atom interaction. In this case, the generator corresponds to a so-called tower feature map. The ratio of the drive amplitude to interaction the strength is^ / ^ ≫ 1. It can be seen that in this case, most spectral gaps are contained approximately inthe interval ^ = [0, 8]. The upper bound may be chosen, e.g., by eye, heuristically, or basedon s predefined criterion, e.g., the smallest interval (either with a fixed lower bound, such as 0, or with a variable lower interval) such that a specified fraction of spectral gaps is containedin that interval. Comparing it to the error function plots in Fig.7 it is evident that for ^ = 4 and^ = 8 error values for (most of) the interval ^ do not exceed ∼5 · 10−2, even for the largestspectral gaps. On the other hand, for ^ = 2, the error starts exceeding 5 · 10−2 already formuch smaller spectral gaps ∼2.5. This translates in poorer estimation of derivative whenusing ^ = 2 (dashed line) which can be confirmed by the plots in Fig. 9A.Fig. 9A shows a computed derivative of a quantum circuit implemented on a(simulated) 8-qubit system with a 2×4 qubit layout, for the generator with the spectrum shownin Fig.8A. The parameters of generator ^^ were selected such that the external driveamplitude Ω is much smaller than the characteristic atom interaction strength ^, i.e., ^ / Ω ≫ 1.This way, the dynamics of the system is mainly determined by the always-on interactionbetween the atoms (qubits). Results are shown for ^ = 2 (dashed line), ^ = 4 (dash-dottedline) and ^ = 8 (dotted line).Fig. 8B shows a distribution of spectral gaps for the same neutral-atom system andthe same tower feature map generator, but in the case that the drive amplitude Ω is close toor slightly larger than the characteristic interaction strength ^, i.e., ^ / ^~1. It can be seen thatthe distribution of gaps is mostly confined to interval of roughly ^ = [0,4] and comparison witherror functions in Fig. 7 shows that already ^ = 2 is sufficient to approximate the derivativewell since the error is sufficiently small. This is confirmed by derivative plots in Fig. 9B whereestimations with all ^ values practically coincide with the exact derivative. Fig.9B showsresults analogous to those shown in Fig.9A, but for the generator with the spectrum shownin Fig.8B. In Figs.9A and 9B, a function ^(^) was calculated using a 0-state in thecomputational basis of the ^-qubit system. With a random initial state, the estimatedderivatives for different values of ^ are given in Figs.9C and 9D.In particular, Fig.9C shows a comparison of derivatives d^ / d^ computed for differentvalues ^ with random initial state versus the exact derivative, wherein the ratio of the driveamplitude to interaction the strength is ^ / ^ ≫ 1. Fig.9D shows a comparison of derivativesd^ / d^ computed for different values ^ with random initial state versus the exact derivative,wherein the ratio of the drive amplitude to interaction the strength is ^ / ^~1. These figures show that in case of a small drive amplitude, the estimation of thederivative is much worse for ^ = 2 (dashed line), relatively acceptable for ^ = 4 (dash-dotted line), and very good for ^ = line). The deteriorated performance ^ = 2 and ^ = 4estimations can be attributed to the fact that since these ^ values produce quite narrowsmall-error intervals, they are not sufficient to cover the tail of distribution in Fig.8A. Withrandom initial state ^^ values from eq. (8) become different and the tail of spectral gapdistribution might become important because modified ^^can increase / decrease the importance of spectral gaps in distribution ^(Δ), thus effectively changing it. A similar (butmuch smaller) effect can be observed in Fig. 9D where the drive amplitude is large (similar inamplitude to the interaction strength). In this case, the spectral gap distribution is narrower,but still, for ^ = 2, the estimated derivative does not coincide with the exact one perfectly andthere is a slight deviation. A much used feature map in quantum machine learning (QML) is the so-called tower feature map, which is used as a part of the quantum circuit to encode a classical (input) variable. Its generator on neutral-atom devices resembles the one in eq.(32), however now the rotation terms ^^^^have distinct prefactors ^^: ^ ^^^ Fig. 10 shows the qubits as was used for Fig.8 and 9, where for qubits with indices ^ ≤ 4 the prefactor is ^^ =1, and for the rest ^ = 2.81. This distribution sh ^ows several prominent clusters of gap starting with 0.4, 5.9, etc. Given that the majority of spectral gaps are found around these peak values and having in mind the structure of the error function in eq. (25), it isevident that using the peak values in algorithms according to embodiments is beneficial sincethe error terms tend to get smaller. Derivative estimations with gaps picked from the prominent peaks of distribution areshown in Fig.11. From these results it can be seen that with ^ = 3 when we select gapscorresponding to the peak values of first three clusters, the estimated derivative fits the exactcurve very closely. In comparison, also with ^ = 3 but selecting equidistant gaps we see thatestimated derivative is not that accurate. With ^ = 2 although gaps are selectedcorresponding to peak values of first two clusters, the resulting derivative is quite far from the exact result. This demonstrates that when having full information about the spectral gaps of some generator it is beneficial to select the gaps according to the most common ones, thus leading to better accuracy of the derivative in comparison with naive equidistant gap sampling.The distribution of spectral gaps shown in Fig. 10A features several clearlyidentifiable clusters of spectral gaps. However, as is shown in Fig.11, not all clusters seem to have similar effect on the accuracy of derivative estimation. The following shows how the clusters may be determined that are the most important to and have the biggest influence on the derivative estimation accuracy. To this end, the distribution of spectral gaps is analysedwhen the interaction term ^ in eq. (33) has been set to zero. Such distribution is shown inFig.10B. Fig. 10B shows that the spectral gap distribution reduces to sharp delta function-likepeaks. The positions of these peaks approximately correspond to the (broadened) cluster peaks in Fig.10A. This observation suggests that information obtained from the rotation part of the generator (without interaction) is still useful for investigating the distribution of the interacting case. However, in the non-interacting case, there are also quite a few sharp peaks in the distribution and it is not immediately clear which peaks (spectral gaps) are the most important. An initial assumption can be to select gaps corresponding to the highest peaks, i.e., in Fig.10B: ∼2, ∼3.5 and ∼5.6. In the current example, it is possible to calculate analytically the expression of thefunction ^(^) corresponding to generator in eq. (33) in the non-interacting case. With thegiven parameters values, the calculated function is a sum of cosine functions with periodsdefined by spectral gaps Δ^ = 2 and Δ^ = 5.63. Thus, only these two gaps are relevant in thenon-interacting case. Since with increasing interaction spectral gap peaks begin to broadenand shift to the right, the peaks in Fig.10A corresponding to gap values 2 and 3.5 merge intoa single peak with gap value ∼3, and the peak corresponding to value 5.63 shifts to around 6. Thus, choosing gap values around 3 and 6 should give an accurate estimation of the derivative, as they correspond to the most important gaps of the spectral distribution. Results in Fig.11B confirm this assumption. The derivative estimated with ^ = 2optimal gap values (dashed line) fits the exact result very closely. The result with twoequidistant gap values (dash-dotted line) is much worse. To achieve a similar result as the 2optimal gap values with equidistant gap sampling, ^ = 4 gap values are needed (dotted line).Thus, using information of the non-interacting generator the most important gaps in the interacting case may be selected. This leads to conclusion that even for a large system with a generator that incorporates an interaction term, it is not always necessary to calculate thefull spectral distribution. Instead, an informed selection of gap values can be made based onthe structure of (only) the non-interacting (here, rotational) part of the generator. This task is usually much easier since the non-interacting generator is often structured and the spectrum can be inferred analytically. To summarize, increasing the number ^ of equations in the linear system (11) resultsin a more accurate estimation of the derivative by reducing the bias error. Compared to thefull GPSR method described by Kyriienko et al., the number ^ of relevant gaps for derivativeestimation can be decreased exponentially, especially when considering larger systems witha number of qubits ^ > 10. The behaviour of the error function ^^ defined in eq. (25)provides a theoretical basis for this drastic reduction of relevant spectral gaps. Moreover, there is no need to consider the exact distribution of spectral gaps since uniform samplingwith, e.g., step 1 from interval [0, ^] works well and reduces the error of gap estimation tonearly 0 in a sufficiently large region. Thus, the embodiments described herein allow tocompute (estimate) derivatives of functions obtained with systems containing large numberof qubits for which calculating all the gaps becomes unfeasible. In the second case when the system size allows calculation of the full spectrum the distribution of gaps can be used to select gap values that correspond to the peak values of the most prominent clusters of gap values. Finally, spectral information obtained from a non-interacting part of the generator may be used to make a smart selection of gaps around the values that are shown to be the most important in the non-interacting case. Fig. 12A–D and Fig.13 illustrate a comparison between estimation of a functionderivative according to an embodiment on the one hand and a Finite Difference (FD) method on the other hand. For this comparison, as an example, a 2-qubit system is considered with a neutral-atom Hamiltonian as a generator that is given by: ^^ ^Here, ^ represents the laser intensity driving the 2-qubit system. The evolution operator of this system can be written in the following way: ^ ^ ^where the differentiation parameter ^ = being the duration when laser pulse with amplitude Ω is applied to the system.Using eq. (7), the analytical expression of function ^(^) can be written explicitly usingthe evolution operator in eq. (35) and the total magnetization cost operator ^^ = ^^^ ⊗ ^^^ +^^^ 2^ and also 4 unique positive spectral ^ In this comparison only aincluded in the linear system of eq. (9), so ^ = 1) is used for the method according to theembodiment. Using the gap value selection procedure described above, the optimal gap thatminimizes the higher order errors is Δ^for the present case. Thus, the corresponding expression of the derivative in this case is given by: ^^ + ^ − ^ ^ −Since the analytical , respectively, the bias ^^,^^(labelled here with “sp” for “spectral”; the plots use the label “aGPSR”) of the derivative estimation can be explicitly computed as the squared difference between the exact derivative and its approximation. This function depends on both thedifferentiation parameter ^ and the parameter shift ^:^ ^The derivative independent random variables ^(^ + ^) and ^(^ − ^), and therefore its variance ^^^,^^is given by: ^^^ ^ ^^ − ^^ of where the variance ^^ ^^^ = ^^ ^ ^^ − ^^For the FD method is given by: ^− ^ ^ −^ and the variance is ^^ ^^ + + ^^ ^ − Based on these resulting from application of the method according to an embodiment (labelled as aGPSR inthe plots) and using a finite difference method (labelled with FD in the plots) can becompared. Fig. 12A–D were generated using generator parameter values Ω = 10 and ^ = 6,number of measurements ^^^^^^ = 10, and differentiable parameter ^ = 2. The comparison isgiven in Fig.12A. This figure shows that for finite differencing, the bias grows more rapidlywith increasing shift ^ compared to method according to an embodiment. In both cases, thevariance sharply grows when shift approaches 0. The figure also shows that the total squared errors, defined as =^ ^ and =^ ^are convex functions weights may be associated with the bias and variance. The respective minima ^^,^^and ^^,^^are points that give an optimal trade-off between bias and variance. For current valuesof Ω and ^, the minimal total error for the finite difference method is larger than that of themethod according to an embodiment. The comparison above is limited to methods that use a single (symmetric) parameter shift; however, with such a simple model, the bias error might be too large for certain applications. The bias error may be reduced by increasing the number of parameter shifts.Fig. 12B–D show a comparison using 2, 3 and 4 parameter shifts, respectively, with therespective minima indicated with ^^,^^and ^^,^^for the finite difference and spectral case,respective, with ^ = 2, 3 ,4. In the current example, the spectral method with 4 parametershifts is equivalent to the (full) GPSR as described by Kyriienko et al., as the system underconsideration has only 4 distinct (or unique) spectral gaps. Consequently, ^^^,^^^^^ = 0 (andnot visible in the plot), = The for the K-shift method according to an embodiment equals, using eq. (14):^^ et^^^^ d(^×)^^ =^^^(^ + ^ ) − ^(^ − )et(^)^ ^ ^(48) d^×^^ ^,^^^^ where the that the ratio of determinants linearly depends ondifferences between shifted function values ^(^ + ^^) − ^(^ − ^^). Hence, ^^,^ is a rationalfunction comprising various combinations of sine functions found in the ^^×^matrix. It isassumed that the optimal spectral gaps {Δ^} are used, corresponding to the generator fromeq. (34). Using eq. (48), it is straightforward to define the corresponding bias and variancefunctions for the K-shift method to an embodiment :Using the ,shift aGPSR and FD methods can be compared, as shown in Fig.12B–D. For the finite differencing method, equispaced central finite differencing is used as described in B.Fornberg , “Generation of Finite Difference Formulas on Arbitrarily Spaced Grids”,Mathematics of Computation 51:184 (1988), pages 699–706,. Here the values for Ω, ^ and ^are the same as in Fig.12A (the 1-shift case). The figures show that the bias error of themethod according to an embodiment is essentially 0, thus not contributing to the totalsquared error. In all cases the minima of the total error are smaller for the aGPSR method compared to the FD method. This suggests that the shift optimization procedure to achieve lowest possible variance in the aGPSR case helps reduce the total squared error. Investigation of total squared error dependence on shift values reveals that for a particular set of generator parameters this error is minimized and the corresponding minimum in the aGPSR case is smaller than that of FD method. However, for a fuller picture of behaviour of the total error using the aGPSR and FD methods, the dependence of the totalerror on the interaction strength ^ may be explored.Fig. 13 shows the total error ratio between FD and aGPSR. For the single parametershift, using the aGPSR method results in a smaller total error of the derivative estimation when compared to single-shift FD approach for a large range of interaction strength values. The graph shows the total error ratio for parameter values ^ = 2, Ω = 10, and 0 ≤ ^ ≤20. Fig. 14A and 14B show results of a method according to an embodiment for ageneral 8-qubit system. In particular, Fig. 14A shows a distribution of spectral gaps for 8qubit system with an arbitrary qubit arrangement when ^ / Ω~1. Fig.14B shows a comparisonof derivatives d^ / d^ calculated for different values ^ with 0 initial state vs exact derivative.The ratio of the drive amplitude to interaction strength is ^ / Ω~1. Calculations of derivatives for different values K when the arrangement of qubits is arbitrary is presented in Fig.14B. The corresponding distribution of spectral gaps is shown inFig. 14A. Good fit of estimated derivative in Fig.14B shows that the underlying arrangementof qubits does not affect the applicability of aGPSR – the increased number K leads to betterderivative estimates. Fig. 21A and 21B are graphs showing a required number of gaps and a relative error,respectively, as a function of system size. The graphs are based on simulations, whichsimulations are performed for a neutral-atom type Hamiltonian. The atoms are arrangedlinearly, on an evenly spaced register. The distance between the atoms is set at 7 µm, withRydberg level 60; this essentially means that the interaction between the atoms is strong, but not so strong that it dominates everything else. In particular, Fig.21A is a graph in which the horizontal axis represents the size ofthe system as a number of qubits, and the vertical axis represents the number of spectralgaps. The black line 2102 represents the total number of gaps in the system, which would allbe required if the known GPSR method is used to calculate the derivative. The grey line 2104 represents the number of gaps required to calculate the derivative with a relative error of 1% using a method according to an embodiment (aGPSR). Note that the vertical axis islogarithmic; the number of gaps required by GPSR grows significantly faster (approximatelyexponentially with the system size) than the number of gaps required by the methods asdescribed herein (aGPSR), which is substantially constant and mostly independent of systemsize. Fig. 21B is a graph in which the horizontal axis represents the size of the system as anumber of qubits, and the vertical axis represents the relative error for different numbers ^ ofspectral gaps used in a method according to an embodiment. Results are shown for 4 (black curve 2112), 6 (dark grey curve 2114), 8 (medium grey curve 2116), and 16 (light grey curve2118) spectral gaps. If more gaps are used, the error decreases; however, the error is mostlyindependent of the system size. Hence, the number of spectral gaps used can be based on a desired accuracy. In this context it is noted that due to the quantum-mechanic processes used in the calculations, and practical limitations on the number of shots, each measurement already has an intrinsic uncertainty. Thus, the number of gaps can be selected based on an inherent uncertainty of the computations. Fig. 22A and 22B are graphs that are analogous to Fig.21A and 21B, respectively.However, instead of a Rydberg interaction, the simulation use a so-called “XY” interaction. Itcan be seen that the results are qualitatively very similar, showing that the same considerations apply for other kinds of interactions, as present in, e.g., systems withsuperconducting qubits, photonic qubits, spin-based qubits, et cetera (results not shown).Fig.15A is a hardware-level schematic of the actions effectuating the logicaloperations shown in circuit diagrams such as Fig.1, 6, and 16–19. Unitary operators, e.g.those to encode the quantum (kernel) feature map and derivatives thereof, can bedecomposed into a sequence of logical gate operations. These logical gate operations are transformations in a quantum Hilbert space over the qubits encoding an input parameter. In order to transform the internal states of these qubits, a classical control stack is used to send pulse information to a pulse controller that affects one or more qubits. The controller may send a sequence of such pulses in time and for each qubit independently. An initializationpulse is used to initialize the qubits into the |0^ state 1502. Then, for example a series ofsingle-qubit pulses is sent to the qubit array in 1504, which may apply a single-layer feature map. Additionally, two-qubit pulse sequences can be used to effectively entangle multiple qubits with a feature map 1506. The duration, type, strength and shape of these pulsesdetermine the effectuated quantum logical operations. 1508 indicates a ‘break’ in thedepicted timeline, which means the sequence of gates may be repeated in a similar fashion in the direction of the time axis 1512. At the end of the pulse sequences, one or more of the qubits are measured 1510. Fig.15B is a hardware-level schematic of the actions effectuating the logicaloperations shown in circuit diagrams such as Fig. 1, 6, and 16–19, specified for aphotonic / optical quantum processor. Unitary operators, e.g. those used to encode thequantum (kernel) feature map and derivatives thereof, can be decomposed into a sequenceof optical gate operations. These optical gate operations are transformations in the quantum Hilbert space over the optical modes. In order to transform the internal states of these modes, a classical control stack is used to send pulse information to a pulse controller that affects one or more modes. The controller may formulate the programmable unitary transformations in a parametrised way. Initially the modes 1514 are all in the vacuum state 1516, which are then squeezed toproduce single-mode squeezed vacuum states 1518. The duration, type, strength and shape of controlled-optical gate transformations determine the effectuated quantum logical operations 1520. At the end of the optical paths, one or more modes are measured with photon-number resolving, Fock basis measurement 1522, tomography or threshold detectors. Fig.15C is a hardware-level schematic of the actions effectuating the logicaloperations shown in circuit diagrams such as Fig. 1, 6, and 16–19, specified for a Gaussianboson sampling device. Unitary operators, e.g. those used to encode the quantum (kernel)feature map and derivatives thereof, can be decomposed into a sequence of optical gate operations. These optical gate operations are transformations in the quantum Hilbert space over the optical modes. In order to transform the internal states of these modes, a classical control stack is used to send information to optical switches and delay lines. The controller may formulate the programmable unitary transformations in a parametrised way. Initially the modes 1526 are all in a weak coherent state, which is mostly a vacuumstate with a chance of one or two photons and negligibly so for higher counts. Subsequently,the photons travel through optical waveguides 1528 through delay lines 1530 and two-modecouplers 1532 which can be tuned with a classical control stack, and which determines theeffectuated quantum logical operations. At the end of the optical paths, one or more modes are measured with photon- number resolving 1534, or threshold detectors. Fig.15D is a hardware-level schematic of the actions effectuating the logicaloperations shown in circuit diagrams such as Fig. 1, 6, and 16–19, which may be executedon a neutral-atom-based quantum computer. On this type of hardware, unitary operators, e.g. those used to encode the quantum feature map and derivatives thereof, can be decomposed in two different kinds of operations: digital or analog. Both of these operations are transformations in the quantum Hilbert space over atomic states. Schematic (a) of Fig.15D depicts a digital quantum circuit 1538, wherein local laserpulses may be used to individually address neutral atoms to effectuate transitions between atomic states which effectively implement sets of standardized or ‘digital’ rotations on computational states. These digital gates may include any single-qubit rotations, and a controlled-pauli-Z operation with arbitrary number of control qubits. Additionally, such digital gate operations may also include 2-qubit operations. Schematic (b) of Fig.15D depicts an analog mode 1546 of operation, wherein aglobal laser light pulse may be applied to groups of, or all, atoms at the same time, with certain properties like detuning, Rabi frequencies and Rydberg interactions to cause multi-qubit entanglement thereby effectively driving the evolution of a Hamiltonian 1544 of theatomic array in an analog way. The combined quantum wavefunction evolves according toSchrödinger’s equation, and particular, unitary operators ^^ = ℋ^ ^, where ℋ^ denotes the Hamiltonian and ^ the time, can be designed by pulse- parametrised coefficientsof the Hamiltonian in time. This way, a parametric analog unitary block can be applied, which entangles the atoms and can act as a variational ansatz, or a feature map, or other entanglement operation. The digital and analog modes can be combined or alternated, to yield a combinationof the effects of each. Schematic (c) of Fig.15D depicts an example of such digital-analogquantum circuit, including blocks 15461-3of digital qubit operations (single or multi-qubit) and analog blocks 15481-3. It can been proven that any computation can be decomposed into a finite set of digital gates, including always at least one multi-qubit digital gate (universality of digital gate sets). This includes being able to simulate general analog Hamiltonian evolutions, by using Trotterization or other simulation methods. However, the cost of Trotterization is expensive, and decomposing multi-qubit Hamiltonian evolution into digital gates is costly in terms of number of operations needed. Digital-analog circuits define circuits which are decomposed into both explicitly-digital and explicitly-analog operations. While under the hood, both are implemented as evolutions over controlled system Hamiltonians, the digital ones form a small set of pre-compiled operations, typically but not exclusively on single-qubits, while analog ones are used to evolve the system over its natural Hamiltonian, for example in order to achieve complex entangling dynamics. It can be shown that complex multi-qubit analog operations can be reproduced / simulated only with a relatively large number of digital gates, thus posing an advantage for devices that achieve good control of both digital and analog operations, such as neutral atom quantum computer. Entanglement can spread more quickly in terms of wall- clock runtime of a single analog block compared to a sequence of digital gates, especially when considering also the finite connectivity of purely digital devices. Further, digital-analog quantum circuits for a neutral quantum processor that are based on Rydberg type of Hamiltonians can be differentiated analytically so that they can be used in variational and / or quantum machine learning schemes, including the differential quantum circuit (DQC) schemes as described in this application. In order to transform the internal states of these modes, a classical control stack is used to send information to optical components and lasers. The controller may formulate the programmable unitary transformations in a parametrised way. At the end of the unitary transformations, the states of one or more atoms may be read out by applying measurement laser pulses, and then observing the brightness using a camera to spot which atomic qubit is turned ‘on’ or ‘off’, 1 or 0. This bit information across the array is then processed further according to the embodiments. Fig.15E is a hardware-level schematic of the actions effectuating the logicaloperations shown in Fig.1, 6, and 16–19, specified for a photonic / optical quantumprocessor. The quantum model can be decomposed into a sequence of optical gate operations. These optical gate operations are transformations in the quantum Hilbert space of the photons. In order to transform the internal states of these photons, a classical control stack is used to send information to a universal multiport interferometer. The controller may formulate the programmable unitary transformations in a parameterized way.Initially the photons 1555 are in Fock states, weak coherent states or coherent states.The duration, type, strength and shape of controlled-optical gate transformations determine the effectuated quantum logical operations 1556. At the end of the optical paths, the modes are measured with photon-number resolving, Fock basis measurement 1557, tomography or threshold detectors. The hardware implementation shown in Fig.15A–15E may also be used toimplement the quantum circuits shown in Fig.16–19, as applicable.Fig.16A depicts a system describing a Digital-Analog implementation of a quantumalgorithm. In particular, Fig.16A represents a Digital-Analog implementation where the qubitinteractions can be switched on and off and there exists single-qubit addressability for the application of single-qubit gates, a process described in 1600. In order to transform the internal states of these qubits, a classical control stack may be used to send information to a controller that affects one or more qubits. The controller may send multiple information carriers (e.g., pulses) in time and for each qubit independently. An initialization protocol isused to initialize the qubits, in this example into the |0^ state 1602. The current exampleshows only four qubits, but the same principles may be extended in a straightforward manner to systems with more (or less) than four qubits. Then, a quantum feature map is applied 1604 to encode quantum information, e.g., an input variable, into a Hilbert space associated with the quantum processor. Following application of the feature map, a variational Ansatz is applied, implemented as a variational quantum circuit 1606. In the current example, the variational Ansatzcomprises three single-qubit gates ^^ – ^^ – ^^ with different rotational angles ^^. Thesesingle-qubit gates are typically sets of standardized or ‘digital’ rotations on computational states applied to different qubits with different rotation angles / parameters. These digital gates include any single-qubit rotations according to the ^^argument of the rotation. A single rotation / single-qubit gate is not sufficient to perform arbitrary rotation regardless of the parameter / angle since it can only rotate over a single axis. In order to accomplish a general rotation block, three gates in series with different rotation angles / parameters can be used, in this example represented by a block of single-qubit gates 1606. Then, the entanglement in this Digital-Analog approach is generated by a wavefunction evolution 1608, described by the block e^^ ℋ ^. During this evolution, the qubits are interacting amongst themselves, letting the system evolve for a specified amount of time. This process produces the necessary entanglement in the system. The combined quantum wavefunction evolves according to Schrödinger’s equation, and particular unitary operators^^ = e^^ ℋ ^ in which ℋ is the Hamiltonian of the system (for example, for neutral atoms the Hamiltonian that governs the system is ℋ =^^ ^^ ∗ ^^ with ^ being the state occupancy of the Rydberg atoms), and ^ is the evolution time. In this way, a parametricanalog unitary block can be applied, which entangles the atoms and can act as a variational Ansatz. After the evolution of the wavefunction 1608, another set of single-qubit gates may beapplied, similar to the process described in the block of gates ^^ – ^^ – ^^ 1606.Then, the wavefunction may be evolved once more, and finally a measurement in the computational basis occurs as described in 1610. Optionally, additional steps may be included, represented by the ellipses, e.g., before,after, or in between the blocks shown. For example, a different initial state than |0^ may beprepared prior to application of the feature map. Furthermore, the blocks can be executed ina different order, for instance, in some embodiments, block 1608 might precede block 1606.One or more of the blocks 1604–1608, or variations thereof, may be repeated one or moretimes prior to the measurement 1610. Fig. 16B represents a typical quantum feature map 1612 that starts by preparing thequantum register to the |0^ quantum state 1602 and then applying single-qubit rotations in allof the qubits with (typically) different rotation angles 1604. In this example, single qubitrotations 1604 (here chosen as ^^(^^(^)) act on each qubit individually and are parametrizedby a variable ^. For a non-linear feature map encoding such as the Chebyshev feature map,the variable ^ may be used as an angle of rotation.The application 1604 of one rotation operation to each qubit in a register may bereferred to as a single layer of rotation operations. This type of encoding of classical information into quantum states is known as angle encoding, which means that a data vectoris represented by the angles of a quantum state. Angle encoding can be found in manyQuantum Machine Learning algorithms, with the main advantage being that it only requires^ = log ^^ qubits to encode a dataset of ^ inputs with ^ features each. Thus, consideringan algorithm that is polynomial in ^, it has a logarithmic runtime dependency on the data size. Fig.16C illustrates different pulses (potentially differing in pulse duration and / or amplitude / height) that can be sent to the quantum circuit (different qubits) to encode theclassical information. The angle of the rotation that is described by the ^^(^) dictates theduration and geometry of the pulse. Fig. 17A presents a quantum circuit based on a Digital-Analog implementation 1700that has no single-qubit addressability and uses a wavefunction evolution that can be turnedon and off for the entanglement generation. Fig.17A differs from Fig. 16A in that it shows atechnique known as data re-uploading that alters the quantum feature map 1704 and 1710.Loading classical data into a quantum system is a non-trivial task, and uploading large amounts of data makes it an even more critical issue. However, when the case of big-data is not considered, techniques like data re-uploading can be utilized. As the name suggests, this entails uploading the classical data multiple times, and encode them via rotational angle / parameter encoding into single-qubit rotations throughout the quantum circuit. Such a data re-uploading technique may start by setting the quantum register to an initial state, e.g., the all-zero state 1702, followed by angle encoding that is based on a feature map 1704, in this case a tower feature map. Following that, as already explained above with reference to 1606, three single-qubit gates in series are used to perform a general rotation block that acts as the Ansatz 1706. Inthis example the encoding is done through the ^^ variable of the ^^ – ^^ – ^^ rotations; however, any type of combined rotations of the (^^, ^^, ^^) could be applied here. Since there is no single-qubit addressability in this example, each ^^and ^^gate is applied at all qubits simultaneously with the same rotation angle ^^. Subsequently, the wavefunction is allowed to evolve for a specified amount of time 1708. Then the classical data are re-uploaded through angle encoding into single-qubit rotations 1710. However, each time the (same) data are uploaded, the information isencoded using different angles than in the previous data uploading steps. For example, theamount of rotation may be doubled in each data (re-)uploading step 1704,1710 (i.e., the rotational angle of the might be increasing as 1–2–4, etc.). Once again, the rotations ^^ – ^^ – ^^ are applied to all qubits simultaneously withthe same rotation angle for each single qubit gate, but different than 1706. The classical information can be repeatedly encoded into the quantum feature map through the data re- uploading technique, while changing the angle of the single qubit rotational angle (1704, 1710, etc.) every time after the wavefunction evolution. The resulting quantum feature map can become more expressive as a tower feature map by doing this process of serial data re- uploading with different levels of rotational angles every time the re-uploading is done. Fig. 17B illustrates a different (and at current technology levels more realistic) form ofa Digital-Analog implementation. In this implementation, the system has single-qubitaddressability, but the wavefunction evolution cannot be turned off, thus staying always on.This system requires less control over the qubit interactions (since they are always on) and only requires single-qubit addressability. The single-qubit gates are applied while the wavefunction evolution (qubit interaction)is also happening. In the shown example, the quantum circuit 1718 starts by encodinginformation in the |0^ quantum state 1702, followed by a quantum feature map 1704(analogous to 1604). Then, a variational Ansatz is applied, comprising both wavefunction evolution e^^ ℋ ^1708 and simultaneous application of a block of single-qubit gates 1716 (similar to block 1606). After multiple repetitions of block 1716, each time with different rotational angles atthe single-qubit gates, the circuit ends with a measurement in the computational basis 1714.Fig.18A illustrates a quantum circuit 1800 that implements a Digital-Analog implementation with single-qubit addressability, wherein the wavefunction evolution (e.g., atom interactions) can be turned on and off. This quantum circuit is equivalent to thequantum circuit in Fig. 17A, with the time of the wavefunction evolution 1808 being ^ = ^ (indimensionless coordinates), while all other operations are the same as the ones described in1700. The same technique of data re-uploading with a tower feature map 1804, 1810 and thesame Ansatz comprising single-qubit gates 1806,1812 followed by wavefunction evolution1808 are used (all similar to 1700).The reason that such a specific time period for the wavefunction evolution 1808 isselected, is due to the fact that such evolution results in the application of Controlled-Z (CZ) gates between pairs of qubits. The benefit of describing the wavefunction evolution as a set of CZ gates lies in the fact that this Digital-Analog implementation of a quantum circuit resembles the structure of the equivalent Digital quantum circuit, which instead of theoperation 1808 performs a set of CNOT gates between pairs of qubits included in thequantum circuit. This allows more straightforward comparisons to be drawn between the Digital and Analog implementation of such feature maps, and to evaluate their relative performance. Fig.18B illustrates a quantum circuit that implements a Digital-Analog implementation with single-qubit addressability 1818, similar to 1800, but with the entanglement generated via CZ gates between the qubits, resembling the atom interaction / entanglement generation. Fig. 19 presents a workflow for generalized circuit differentiation. In QuantumMachine Learning (QML) algorithms, various quantum feature maps and Parametric Quantum Circuits (PQCs) are typically utilised to express the circuit required to solve a particular problem. These quantum circuits contain tuneable parameters that use a classicaloptimizer to find the optimal values of the parameters that exist for that formulation. Thedifferentiation required for most parameter optimization methods can be achieved in various ways, with the most common one being the Parameter Shift Rule (PSR). One way to calculate analytic derivatives of quantum circuits is by measuring overlaps between quantum states. Analytic derivatives can be used for differentiatingunitaries like the one presented in 1904, e.g., ^^ = e^^ ^ ^^ / ^ generated by arbitrary Hermitiangenerator ^^. However, in order to have a less resource-intensive method for differentiation, the parameter shift rule (PSR) was proposed. The PSR algorithm can provide analytic gradient estimations through measurement of the expectation values, with gate parameters being shifted to different values. The PSR algorithm is much used to perform differentiation for QML algorithms; however, it is only valid for a specific type of generators, viz., generators that are involutory or idempotent. This is because the full analytic derivative only requires 2 measurements of expectation values (2 unique eigenvalues). Although this simplifies the differentiation protocol, it also restricts it to only work for a certain type of generators. Therefore, as an improvements and generalizations of the PSR algorithm, in 1900 ascheme for a Generalized Parametric Shift Rule (GPSR) algorithm is shown. Such approach allows for differentiating generic quantum circuits with unitaries generated by operators with a rich spectrum (i.e., not limited to idempotent generators). A non-exhaustive list of situations when such generators arise in the implementations shown in Fig.15–18: when the parameter controls a multi-qubit gate, when the single-qubit is not perfectly isolated from other qubits due to interaction or when the parameter can only simultaneously control several single qubit gates. The GPSR algorithm is based on spectral decomposition, and showcases the role of the eigenvalue differences (spectral gaps) during differentiation for the generator spectrum.Such approach works for multiple non-degenerate gaps, in contrast to the PSR that is onlyvalid for involutory and idempotent operators with a single unique spectral gap. A workflow for generalized circuit differentiation can be described as follows: create aquantum circuit that encodes the function ^ similar to 1900 including single- and multi-qubitgates 1902. Pick a unitary operator such as ^^(^) = exp^−i ^ ^^ / 2^ 1904 that is parametrizedby some tuneable parameter ^, and study the spectrum of its generator ^^. Then, using unique and positive spectral gaps, a system of equations that includes the spectral information can be created, alongside the calculated parameter shifts, and the measured function expectation values as shifted parameters. The solution of such system can providethe analytical derivative for any general generator ^^ and thus perform GPSR. GPSR isdescribed in more detail in Kyriienko et al., ‘Generalized quantum circuit differentiation rules’,Phys. Rev. A 104 (052417), which is hereby incorporated by reference.The generalized parameter shift rule provides a method to calculate the exact derivative with respect to parameters represented by an arbitrary generator with a rich spectrum of eigenvalues while remaining robust against noise. However, the generator of aquantum circuit with ^ qubits will have 2^ Eigenvalues, which may all be distinct, whichcorresponds to making measurements at up to 2^(2^ − 1) different values of the parameterto estimate the derivative at any point. This quickly becomes prohibitively expensive. Theembodiments described herein provide methods for computing (approximating) the derivativeof a quantum circuit with a limited shot budget, in a system with a non-trivial spectrum of Eigenvalues (i.e., more than two distinct Eigenvalues, but typically hundreds or even thousands). The workflow and the area of applicability remains the same as 1900, thedifference being in the number and selection of the gap values that are used in thecalculation. Fig. 20A depicts a system for determining a solution for an computational problemusing a hybrid quantum computer system according to an embodiment. The system 2002may include a quantum computer system 2004 comprising one or more quantum processors2008, e.g. a gate-based qubit quantum processor, and a controller system 2010 comprisinginput output (I / O) devices which form an interface between the quantum processor and the outside world, e.g., the one or more classical processors of a classical computer 2006. For example, the controller system may include a system for generating control signals for controlling the quantum processing elements. The control signals may include for example a sequence of pulses, e.g. microwave pulses, voltage pulses and / or optical pulses, which are used to manipulate qubits. Further, the controller may include output device, e.g. readoutcircuits, for readout of the qubits and the control signals for readout of the quantumprocessing elements, e.g. a readout pulse for reading a qubit. In some embodiments, at least a part such readout circuit may be located or integrated with the chip that includes the qubits. The system may further comprise a (purely classical information) input 2012 and an(purely classical information) output 2014. The data processor systems may be configured to solve an optimization problem using the quantum computer. Input data may include information about the optimisation problem one wishes to solve. This information may include training data such as input-output pairs, a differential equation, boundary conditions, initial values, regularization values, etc. The input data may be used by the system to construct quantum circuits, in particular quantum feature maps and / or parametrized quantum circuits,and to classically calculate values, e.g. sequences of pulses, which may be used to initializeand control qubit operations according to the quantum circuit. To that end, the classical computer may include a quantum circuit generator 2007. The input data may be used by the system to classically calculate values, e.g. parameter settings, which may be used to initialize the quantum circuit that is implemented on the quantum processor. Similarly, output data may include loss function values, sampling results, correlator operator expectation values, optimisation convergence results, optimized quantum circuit parameters and hyperparameters, and other classical data. Each of the one or more quantum processors may comprise a set of controllable quantum processing elements, e.g. a set of controllable two-level systems referred to asqubits. The two levels are |0^ and |1^ and the wave function of a ^-qubit quantum processormay be regarded as a complex-valued superposition of 2^of these (distinct) basis states. The embodiments in this application however are not limited to qubits but may include any multi-level quantum processing elements, e.g. qutrits, that is suitable for performing quantum computation. Examples of such quantum processors include noisy intermediate-scale quantum (NISQ) computing devices and fault tolerant quantum computing (FTQC) devices. The quantum processor may be configured to execute a quantum algorithm in accordance with the gate operations of a quantum circuit. The quantum processor may be implemented as a gate-based qubit quantum device, which allows initialization of the qubits into an initial state, interactions between the qubits by sequentially applying quantum gates between different qubits and subsequent measurement of the qubits’ states. To that end, the input devices may be configured to configure the quantum processor in an initial state and to control gates that realize interactions between the qubits. Similarly, the output devices may include readout circuitry for readout of the qubits which may be used to determine a measure of the energy associated with the expectation value of the Hamiltonian of the system taken over the prepared state. Furthermore, the one or more quantum processors may comprise a set of continuous variable systems, such as optical or photonic quantum computers. Furthermore, the quantum processor may comprise a Gaussian Boson Sampling device as described above withreference to Fig. 15C or any of the other hardware described above with reference to Fig.15A–E. In some embodiments, the first data processor system may be implemented as asoftware program for simulating a quantum computer system 2004 comprising a quantumprocessor system 2008. Hence, in that case, the software program may be a classicalsoftware program that runs a classical computer 2006 so that quantum algorithms can bedeveloped, executed and tested on a classical computer without requiring access to a hardware implementation of the quantum processor system. Fig.20B is a block diagram illustrating an exemplary classical data processing system described in this disclosure, for example classical computer 2006. Classical dataprocessing system 2020 may include at least one processor 2022 coupled to memoryelements 2024 through a system bus 2026. As such, the data processing system may storeprogram code within memory elements 2024. Further, processor 2022 may execute theprogram code accessed from memory elements 2024 via system bus 2026. In one aspect,data processing system may be implemented as a computer that is suitable for storing and / or executing program code. It should be appreciated, however, that data processing system2020 may be implemented in the form of any system including a processor and memory thatis capable of performing the functions described within this specification. Memory elements 2024 may include one or more physical memory devices such as,for example, local memory 2028 and one or more bulk storage devices 2030. Local memorymay refer to random access memory or other non-persistent memory device(s) generally used during actual execution of the program code. A bulk storage device may be implemented as a hard drive or other persistent data storage device. The classical data processing system 2020 may also include one or more cache memories (not shown) that provide temporary storage of at least some program code in order to reduce the number oftimes program code must be retrieved from bulk storage device 2030 during execution.Input / output (I / O) devices depicted as key device 2032 and output device 2034optionally can be coupled to the data processing system. Examples of key device may include, but are not limited to, for example, a keyboard, a pointing device such as a mouse,or the like. Examples of output device may include, but are not limited to, for example, amonitor or display, speakers, or the like. Key device and / or output device may be coupled to data processing system either directly or through intervening I / O controllers. A networkadapter 2036 may also be coupled to data processing system to enable it to become coupledto other systems, computer systems, remote network devices, and / or remote storage devices through intervening private or public networks. The network adapter may comprise a data receiver for receiving data that is transmitted by said systems, devices and / or networks to said data and a data transmitter for transmitting data to said systems, devices and / or networks. Operation modems, cable operation modems, and Ethernet cards are examples of different types of network adapter that may be used with classical data processing system 2020. As pictured in FIG. 20B, memory elements 2024 may store an application 2038. Itshould be appreciated that classical data processing system 2020 may further execute anoperating system (not shown) that can facilitate execution of the application. Application, being implemented in the form of executable program code, can be executed by classical data processing system 2020, e.g., by processor 2022. Responsive to executing application, data processing system may be configured to perform one or more operations to be described herein in further detail. In one aspect, for example, classical data processing system 2020 may represent aclient data processing system. In that case, application 2038 may represent a clientapplication that, when executed, configures classical data processing system 2020 toperform the various functions described herein with reference to a “client”. Examples of a client can include, but are not limited to, a personal computer, a portable computer, a mobile phone, or the like. In another aspect, data processing system may represent a server. For example, data processing system may represent an (HTTP) server in which case application 2038, when executed, may configure data processing system to perform (HTTP) server operations. In another aspect, data processing system may represent a module, unit or function as referred to in this specification. The terminology used herein is for the purpose of describing particular embodiments only and is not intended to be limiting of the disclosure. As used herein, the singular forms “a,” “an,” and “the” are intended to include the plural forms as well, unless the context clearly indicates otherwise. It will be further understood that the terms “comprises” and / or “comprising,” when used in this specification, specify the presence of stated features, integers, steps, operations, elements, and / or components, but do not preclude the presence or addition of one or more other features, integers, steps, operations, elements, components, and / or groups thereof. The corresponding structures, materials, acts, and equivalents of all means or step plus function elements in the claims below are intended to include any structure, material, or act for performing the function in combination with other claimed elements as specifically claimed. The description of the present disclosure has been presented for purposes of illustration and description, but is not intended to be exhaustive or limited to the embodiments in the form disclosed. Many modifications and variations will be apparent to those of ordinary skill in the art without departing from the scope and spirit of the invention. The embodiments were chosen and described in order to best explain the principles and the practical application, and to enable others of ordinary skill in the art to understand the various embodiments with various modifications as are suited to the particular use contemplated.

Claims

CLAIMS 1. A method for computing a derivative of a quantum circuit using a hybrid dataprocessing system, the hybrid data processing system comprising a classical computer system and a quantum computer system, the quantum computer system comprising a plurality of quantum elements, the method comprising: determining or receiving, by the classical computer system, a formulation of a parametrized quantum circuit parametrized by a continuous variable in a variable space, the parametrized quantum circuit encoding a given parameter value of the continuous variable in a Hilbert space associated with the plurality of quantum elements; determining or receiving, by the classical computer system, a number of one or moredistinct gap values, the number of distinct gap values being smaller than a number of distinctspectral gaps obtained from an Eigenvalue spectrum of a generator associated with theparametrized quantum circuit; determining or receiving, by the classical computer system, a number of distinctparameter shifts, the number of distinct parameter shifts being equal to the number of one ormore distinct gap values;executing, by the classical computer system, for each of the one or more distinct parameter shifts, the parametrized quantum circuit with the continuous variable shifted basedon the respective parameter shift, the execution comprising:– translating the parametrized quantum circuit into control signals for controlling theplurality of quantum elements and for readout of the plurality of quantum elements to obtain hardware measurement data; –controlling the quantum computer system based on the control signals; and– receiving, in response to the execution of the parametrized quantum circuit, thehardware measurement data; and processing the hardware measurement data to obtain the derivative of the parametrized quantum circuit with respect to the given parameter.

2. The method as claimed in claim 1, further comprising determining or receiving, bythe classical computer system, spectral information describing the Eigenvalue spectrum;wherein the one or more distinct gap values are based on the spectral information.

3. The method as claimed in claim 2, wherein the spectral information comprises atleast one of: –an exact and / or approximate value of a smallest and a largest positive spectralgap, and / or a smallest and a largest Eigenvalue;– a statistical representation of a distribution of positive spectral gaps; and– a list of all positive spectral gaps.

4. The method as claimed in any one of the preceding claims, wherein the one ormore distinct gap values are obtainable by at least one of:– uniform or quasi-uniform sampling in a predetermined range and / or a rangedefined by the spectral information; and –sampling based on the spectral distribution.

5. The method as claimed in claim 2 or 3, wherein the generator is based on a non-interacting part corresponding to a plurality of non-interacting quantum elements and aninteracting part based on an interaction of the plurality of quantum elements, wherein acontribution of the interacting part to the spectral gaps of the generator is relatively small compared to a contribution of the non-interacting part, wherein the spectral information comprises information on the Eigenvalue spectrum of the non-interacting part, and whereinthe one or more distinct gap values are obtainable by at least one of:– uniform or quasi-uniform sampling in a range defined by the spectral gaps of thenon-interacting part; and –sampling based on the distribution of positive spectral gaps of the non-interactingpart.

6. The method as claimed in any one of the preceding claims, wherein the distinctparameter shifts are obtainable by at least one of: –selecting the parameter shifts based on the gap values, preferably selecting theparameter shifts such that the variance of the derivative is minimized; and– selecting, preferably uniformly selecting, the parameter shifts from afrom the ^^ ^^ ^ , ^ ^, more preferably from the7. The method as claimed in any one of the preceding claims, wherein obtaining thederivative of the parametrized quantum circuit with respect to the given parameter comprises: for each of the one or more distinct parameter shifts, determining a difference between the two shifted function values;solving, based on the determined differences, the one or more distinct gap values andthe one or more distinct parameter shifts, a linear system of equations to determine a gap coefficients associated with each spectral gap, preferably solving ^ ^where ^^denotes the the ^thparameter shift ^^, and^ is the number of distinct gap values Δ^ and distinct parameter shifts ^^and computing the derivative based on the one or more distinct spectral gaps and the associated gap coefficients, preferably by computing ^wherein ^^ ^(^) denotes hequantum circuit ^(^) with respectto the parameter ^ parametrizing the parametrized quantum circuit, ^^ denotes the gapcoefficient associated with gap value Δ^, and ^ is the number of distinct gap values Δ^.

8. The method as claimed in any of the preceding claims, wherein the number ofdistinct gap values is smaller than or equal to half the number of distinct spectral gaps, preferably smaller than or equal to a base-2 logarithm of the number of distinct spectral gaps; and / or wherein the plurality of quantum elements comprises at least eight quantum elements and wherein the number of distinct gap values is at most 64, preferably at most 32, more preferably at most 16, even more preferably at most 8.

9. The method as claimed in any one of the preceding claims, wherein the number ofdistinct gap values is based on an relative error associated with the number of distinct gap values.

10. The method as claimed in any one of the preceding claims,wherein executing the parametrized quantum circuit comprises translating theparametrized quantum circuit into a sequence of signals, at least one of the signals beingdepend on the given parameter, and using the sequence of signals to operate the plurality ofquantum elements; and / or, wherein receiving hardware measurement data comprises applying a read-out signal to the plurality of quantum elements and, in response to the read-out signal, measuring quantum hardware measurement data.

11. The method as claimed in any of the preceding claims, wherein the hybridquantum computer includes a gate-based qubit device, a digital / analog quantum device, a neutral-atom-based quantum device, an optical qubit device and / or a gaussian boson sampling device.

12. A hybrid data processing system for computing a derivative of a quantum circuit,the hybrid data processing system comprising a classical computer system and a quantum computer system, the classical computer system comprising at least one processor and a memory storing computer-readable instructions, the memory being communicatively coupled to the at least one processor, the quantum computer system comprising a quantum register with a plurality of quantum elements, the at least one processor being configured to, in response to executing the computer-readable instructions: determine or receive a formulation of a parametrized quantum circuit parametrized bya continuous variable in a variable space, the parametrized quantum circuit encoding a given parameter value of the continuous variable in a Hilbert space associated with the plurality of quantum elements; determine or receive a number of one or more distinct gap values; determine or receive a number of distinct parameter shifts, the number of distinctparameter shifts being equal to the number of one or more distinct gap values; executing, for each of the one or more distinct parameter shifts, the parametrized quantum circuit with the continuous variable increased by the respective parameter shift and with the continuous variable decreased by the respective parameter shift, the execution comprising: –translating the parametrized quantum circuit into control signals for controlling theplurality of quantum elements and for readout of the plurality of quantum elements to obtain hardware measurement data; –controlling the quantum computer system based on the control signals; and– receiving, in response to the execution of the parametrized quantum circuit, thehardware measurement data; and processing the hardware measurement data to obtain the derivative of the parametrized quantum circuit with respect to the given parameter.

13. The system as claimed in claim 12, configured to execute the method as claimedin any one of claims 1–11.

14. A computer program product or suite of computer program products comprisingsoftware code portions which, when run by a hybrid data processing system comprising a classical computer system and a quantum computer system, configure the hybrid dataprocessing system to execute the method steps according to any one of claims 1–11.

15. A non-transient storage medium storing a computer program product according toclaim 14.

Images