Parallel Quantum Computation of Vibrational Dynamics

The vibrational dynamics in a linear triatomic molecule is emulated by a quantum information processing device operating in parallel. The quantum device is an ensemble of semiconducting quantum dot dimers addressed and probed by ultrafast laser pulses in the visible frequency range at room temperature. A realistic assessment of the inherent noise due to the inevitable size dispersion of colloidal quantum dots is taken into account and limits the time available for computation. At the short times considered only the electronic states of the quantum dots respond to the excitation. A model for the electronic states quantum dot (QD) dimers is used which retains the eight lowest bands of excitonic dimer states build on the lowest and first excited states of a single QD. We show how up to 82 = 64 quantum logic variables can be realistically measured and used to process information for this QD dimer electronic level structure. This is achieved by addressing the lowest and second excited electronic states of the QD’s. With a narrower laser bandwidth (= longer pulse) only the lower band of excited states can be coherently addressed enabling 42 = 16 logic variables. Already this is sufficient to emulate both energy transfer between the two oscillators and coherent motions in the vibrating molecule.


INTRODUCTION
We describe the theoretical background for an experimental setup, an ensemble of quantum dot dimers that can and has been realized in the laboratory. We show explicitly how this device is used to emulate the quantum vibrational dynamics of a linear triatomic molecule and discuss possible extensions. In 1985 Deutsch defined a quantum computer as a device that could simulate effectively an arbitrary physical system [1]. Our aim here is much more modest. We seek to describe a device that can be realized with currently available laboratory techniques. Furthermore, the device needs to provide computational answers only for a limited set of variables of the physical system. The computation is realized by mapping of the dynamics of the physical variables of this limited set using a set of observations of the time-evolution of the device. The set of possible observations of the device is the set of our N 2 logic variables. The number N is less than or equal to the number of accessible quantum states of the logic device. N 2 is larger than the number of variables of interest for the physicochemical system that is emulated.
The quantum mechanical device is characterized by a set of observables and that is why we call the logic done by the device "computing by observables" [2,3]. With such a device one can emulate a physical system that requires up to this number for simulating a closed set of physical variables.
Increasing levels of resolution provided by the description of a physicochemical system can be characterized by an increasing numbers of variables. An example used in this paper is that of two coupled harmonic oscillators. For each oscillator there is a denumerable number of bound states. So, in principle, a fully complete characterization of each oscillator can require up to a denumerably infinite number of variables. An opposite extreme is a thermal state of the oscillator. Its complete description requires two physical variables, a normalization of the probabilities and the mean energy. For the coupled oscillators we will describe several physical levels of resolution that are intermediate between these two limits.
The variables that we deal with can be regarded as vectors in a linear space. They can, if so desired, be made orthogonal etc. As vectors, the variables seem to behave like classical variables but they are a quantum mechanical description of the system. In textbooks one is more used to a linear vector space of wave functions, what is sometimes called a Hilbert space. The observables are also defined in a space with a scalar product and a Mathematician may choose to call it a Hilbert space. But if the Hilbert space of the wave functions has N dimensions, the corresponding linear vector space of observables has N 2 dimensions. If {|i〉} is a set of basis vectors in the Hilbert space of wave functions, i 1, 2, .., N, a basis set for the observables can be written as {X ij |i〉〈j }, i, j 1, 2, .., N. The N 2 expectation values of these are our variables. We will refer to the expectation values when i ≠ j as coherences. The N diagonal elements are populations. By taking linear combinations, e.g., |i〉〈j + j〉〈i|, we can arrange all the variables to be Hermitian and all the variables to be real.
In this paper we use two distinct linear vector spaces. One is suitable to describe the dynamics of the computing device and one the dynamics of the physical system we aim to emulate. We do not consider the device as providing an analog for the physical system. To emphasize this point we use a device that operates on the dynamics of purely electronic states pumped and read by ultrafast laser pulses. A schematic description of the device and the physical system is shown in Figure 1. The dynamics in the electronic Hilbert space of the quantum dot dimers is mapped to simulate the time evolving physical variables in the nuclear Hilbert space. We will therefore need a separate discussion of the logic variables provided by the computing device and those required to describe the physical system. Quantum computing algorithms have been developed for problems of chemical interest. This includes not only electronic structure [4][5][6][7][8][9][10][11][12][13] but also quantum enhanced machine learning algorithms [8] and algorithms for dimension reduction including Principal Component Analysis, Canonical Correlation Analysis and other algebraic methods used for dimension reduction such as surprisal analysis [6,[15][16][17] A preliminary report on our quantum computation has been published [18].

THE COMPUTING DEVICE
Our computing device is an array of CdSe quantum dot QD dimers. Addressing the device is by a sequence of laser pulses and read-out is performed simultaneously on many dots [19][20][21][22][23][24] as in 2D electronic spectroscopy [25,26]. Measuring over a classical ensemble of dots importantly means that we directly read the mean values and that there is no interference between measuring different observables that individually do not commute. In the short time interval during which the time and probe are performed the primary source of noise is the variability in the size of the quantum dots. This means that the two dots making the dimer are not quite identical so strictly speaking we have heterodimers. We therefore need to average the read-out over the distribution of dimer sizes. The states of different dimers are at slightly different energies so that the FIGURE 1 | Representation of the computing device (A) and the physical system (B). The four bands of excited electronic states of the device that are used are indicated as well as the ten pairs i, j whose twenty corresponding coherences X ij are used as logic variables. The physical system we emulate is a vibrating triatomic molecule in a non-stationary quantum mechanical state. We will simulate the mean position and momentum of each of the two coupled oscillators as well as the respective widths, σ R and σ P .
Frontiers in Physics | www.frontiersin.org October 2020 | Volume 8 | Article 590699 coherences between different pairs of states beat at slightly different frequencies. This reduces the frequency resolution or, equivalently, the time available for read-out [13] to faster than the dephasing time of the mixture of not fully identical dimers. So we do not read individual quantum states but bands of closely adjacent states belonging to different dimers. In our computations we use an experimentally realizable size of the dots and their size dispersion.
In this paper we employ rather elementary dynamics of the computing device. Using a sequence of laser pulses we bring each dimer to a multicoherent state. From then on the device is unperturbed and it is the oscillations of different coherences that are used to simulate the time dependence of the physical system. The coherence between the ground state and any electronic excited state oscillates quite fast so it dephases faster then the dephasing of the coherences between the electronically excited states. The very fast beating coherence, few fs's or less, between the ground state and an electronic excited state cannot therefore be reliably detected. So very fast oscillations cannot be simulated by our device due to the currently available dispersion in sizes of colloidal quantum dots. It is expected that a lower dispersion will be experimentally possible in the near future.
This paragraph is a sketch of the counting of electronically excited states in a CdSe dot of nanometric mean size between 2 and 5 nm in diameter and of a dimer of two dots drawn from an ensemble of dots of the same mean size. The two dots making up the dimer are of about the same but not identical size. The first electronically excited state of a dot is a band of 12 states, as follows. The spin of the hole is coupled to a p (l 1) orbital of Se to give rise to two spin orbit components, 1/2 and 3/2. The s orbital on Cd contributes only a spin of 1/2. The total orbital angular momentum F of the excitonic 3/2 state has the projections F z ± 2, ± 1, 0 and F z ± 1, 0 for the hole and the electron respectively, eight states all together. There are four states, (F z ± 1, 0 and 0) for the 1/2 excitonic state. The two excitonic bands of states, (of eight and four states, respectively), can be experimentally resolved [23]. There are similarly two resolvable bands of states for the second electronically excited state of CdSe. All together, four bands of electronically excited states per dot. In the dimer, or dots with short ligands, the dots are interacting by Coulomb and exchange coupling. Thereby each one of the four bands of states of the monomer is split into two. In the visible range a dimer has eight resolvable bands of states, four bands that correspond primarily to the lower electronic state and four to the higher excited state of a monomer. Since the dimers are made of two dots that are not quite identical, all eight bands are optically active from the common dimer ground state, see Supplementary Materials, section 1.1 for more details.
Using the counting as above the CdSe dots enable a choice of four or eight bands of states in the visible range. Addressing coherently the states of all the eight bands requires shorter laser pulses. Smaller QD have larger energy difference between their two lower excited electronic states and therefore one needs a larger energy bandwidth of the pulse or a higher carrier frequency to simultaneously address them. On the other hand, for larger dots, the energy difference between the two lowest excited states may become smaller or comparable to strength of the spin-orbit coupling, which leads to a loss of resolution between the different bands [23,24]. Commercially available fast lasers in the visible can easily address coherently all the states of the four lowest bands. The minimal capability of our device is therefore 5 2 logic variables, the five populations of the ground and four excited states and the 5 × 4 coherences. Four (times two, complex values) coherences of the transitions between the ground state and the four excited states and 12 coherences between the four excited states. After the laser pulses are over, each coherence will oscillate with a fixed frequency determined by the differences in energy of the two states. The highest frequencies are for the transitions up from the ground state. These are rather fast. So in this paper we consider the coherences between pairs of excited states.
To describe the dynamics of the device we take it that initially it is in the ground electronic state. We assume that the addressing lasers are weak enough so that only one photon transitions are possible. The lasers are in the visible so that transitions between the excited states are way out of resonance and so are excluded. Indeed and as we shall discuss, the frequencies of the transitions between the excited states are in the range of molecular vibrational frequencies so that we can use the coherences between excited states to emulate the time scales relevant to the physical systems. The range of timescales is determined by the size of the quantum dot dimers as is discussed above and by the coherence width of the excitation lasers. The bandwidth of a 6 fs pulse is about 2,100 cm −1 .
The excitation scheme is that typically used in 2D electronic spectroscopy [25,26]. The first fast laser pulse generates the absorption from the ground state. Next, with some delay, is a second laser pulse. One frequency axis that we will use is the Fourier transform with respect to this time delay. The system is next allowed to evolve for a time interval that is typically denoted T. After the second pulse, the system can be back in the ground state, in one of the excited states, or in a coherence between two excited states. After the time T the third laser pulse stimulates emission that is monitored in time and the second frequency axis is a Fourier transform with respect to this time. So for each value of the time T we generate a 2D frequency map. In the echo phase matching direction of emission, the populations of the excited states appear on the diagonal of the 2D maps and the coherences are on the off-diagonal at a position determined by the excitation frequency of the two excited states they connect. At different values of the time interval T, the intensities at the off-diagonal positions on the map will vary according to the time-evolution of the respective coherences, and this is how we can simulate the time dependence of the physical variables. As emphasized above, due to the size dispersion the coherence beating frequencies will have a finite spread. The ability to resolve coherences between different pairs of states depends on this spread to be limited enough. From the complementary time-dependent point of view, the coherence contributions decay in time as a Gaussian with a width that is governed by the width of the frequency distribution of the coherence [3]. On the plus side, the finite spread of the frequency of every coherence means that the computing device can simulate finite spans of values, roughly the width of the Gaussian, about different frequencies.
Frontiers in Physics | www.frontiersin.org October 2020 | Volume 8 | Article 590699 Deterioration of the signal due to dephasing has been an issue for quantum computing from its early beginning [27,28]. Often it is dephasing due to interaction with the environment. Here too, the limiting factor is dephasing but it is a dephasing due to an inhomogeneous broadening caused by the size dispersion of the QD's. In the short times that we probe, the electronic states of each dimer practically are not yet perturbed. The dephasing is due to the probe averaging over many dimers, dimers with slightly different frequencies [3]. We pay special care to make the dephasing realistic. The mean diameter of the CdSe QD's, 4.4 nm, that we use to compute the spectrum and the size dispersion, 5%, have been realized experimentally [23,29]. We emphasize the considerable improvement in performance that can be achieved by an even modest reduction in the dispersion in sizes. Yet even the already accessible 5% size variability is sufficient to enable a quite powerful device. We also stress that the 2DES measurements that we rely upon are performed at room temperature, in solution [23] or on a solid state device [24].
We can emulate dynamics by measuring 2D coherence maps at a range of values for the time interval T. With current chemical synthetic capabilities the dephasing due to the size dispersion of the dot is the primary limitation on how long T can be. At longer times, say beyond 100 fs, coupling to the phonons also sets in.
For each experiment there will be a sequence of transitions leading to the same coherence between two excited states during the delay time T between the second and the third pulse. For weak pulses so that first order time-dependent theory and the impulsive limit transition applies at each interaction with the pulse, the two coordinates of the 2D maps are the absorption frequency from the ground state (GS) and the emission frequency to GS. When two excited states, i and j, can be reached by one photon transitions from the ground state (GS) there are four positions in the 2D maps, two diagonal ones at (ω 0i , ω 0i ) and (ω 0j, ω 0j ) and two off-diagonal ones (ω 0i , ω 0j ) and (ω 0j , ω 0i ) where ω 0i and ω 0j are the excitation frequencies of states i and j from the GS. In the rephasing phase matching (echo) direction, eight Liouvillian paths contribute to these four positions. They are spelled out explicitly in Supplementary Figure S3 in the notation of ref. [26], see also refs. [2,30].
Four paths fall on the diagonal of the 2D map at the excitation frequencies of the two excited states: one ground state bleaching (GSB) path for which the system is the GS during the time interval T and one stimulated emission paths (SE) for which the system is in excited state i during T at (ω 0i , ω 0i ) and a GSB and a SE path for which the system is in excited state j at (ω 0j , ω 0j ). At the short times below 200 fs considered here, there is no exchange of population between the GS and the two excited states during T and the contribution of these four paths to the 2D response is time-independent. Two time-dependent GSB paths contribute to the two off diagonal positions, (ω 0i , ω 0j ) and (ω 0j , ω 0i ), respectively. Two paths that are the signature of the electronic coherence between the two excited states also contribute at the same off diagonal positions. The two coherence paths are beating during T with the transition frequency ω ij of the coherence between the excited states i and j.
The eight paths, see Supplementary Figure S3, contribute simultaneously to the time evolution of the density matrix and therefore to the 2D map, hence "parallel" in our title. When the system has more than two excited states, each pair of excited states contributes 8 paths to the time evolution of the density matrix. For example, there are 16 coherences between excited states for the five level system (GS + 4 excited states) discussed above, meaning that 16x8 paths contribute in parallel. Of these, only the 16x2 paths that lead to a distinct off diagonal position on the 2D maps are used to map observables of the system that is emulated. There are also 8 paths per pair of excited states in the non rephasing phase matching direction. So in principle, the larger the number of states in the band that can be accessed by one photon transition from the GS, the larger the number of observables that can be emulated.
As explained above, colloidal QD's are synthesized with a finite size dispersion which leads to a distribution of transition frequencies for each coherence, both from the GS, the ω 0i type and between excited states (the ω ij type). At the level of the ensemble, the inherent size dispersion of the QD's leads to a Gaussian distribution in energy of the addresses (ω 0i , ω 0j ) of coherences between excited states and to their Gaussian dephasing along T and therefore limits the number of coherence positions that can be resolved on the map. When controlled and limited to a few percent in diameter, the size dispersion can also be used to advantage. Scanning positions on the 2D maps around the address corresponding to the mean transition frequencies from the GS gives slightly different periods of the coherence along T which provides flexibility in mapping the periods of the observables of the emulated systems. Figure 1 shows a scheme where there are four bands of excited states that are accessed by the laser pulses. The four coherences at a relatively high frequency when the ground state is one of the two states are not detected in 2D spectroscopy. There are the twelve lower frequency coherences that, as mentioned, will be used to emulate the vibrational motions of the physical system. Note that there is a spread in the frequencies and the very lower ones will be used to emulate splittings due to coupling of two vibrational modes. The coherences X ij |i〉〈j| between excited levels i and j of the quantum dot dimer are our logic variables. After the laser pulses are over, the coherences oscillate in time with a frequency that is the difference between the energies of the two states that they connect. We aim to simulate the time dependence of observables of the coupled oscillators system. For example, as shown in Eq. 7, we need the time dependence of 〈a 2 〉 to compute the variance of the bond displacement of the first oscillator.

THE PHYSICAL SYSTEM: LINEAR VIBRATING TRIATOMIC MOLECULE
In this paper we track the electronic dynamics in the QD dimers to emulate the vibrational motion in non-stationary states of a triatomic molecule, Figure 2.
The system we emulate is two coupled harmonic local vibrations, each representing a bond. The coupling is due to the motion of the central atom and it depends on its mass and its displacement from equilibrium [31]. By expanding the potential of the system and keeping only up to quadratic terms in the deviation of the bond displacements from equilibrium, the Hamiltonian for two harmonically coupled harmonic modes, denoted a and b, is, using Z 1, (1) Here we use the standard notation for the creation and annihilation operators for the two local vibrational modes. The vibrational quantum number mismatch is determined by 〈 a † a − b † b〉 and the frequency mismatch is δω ω a − ω b . The harmonic coupling is a b † + a † b and the strength is α/2. Due to this coupling, the creation and annihilation operators of the two oscillators are correlated. We characterize vibrational motion of the local modes by the mean values of the bond distance, momenta and dispersion of both quantities. Computation of the values of the mean bond distance and mean momenta involves description of the evolution of both creation and annihilation operators: 〈 R a 〉 〈 a + a † 〉/ 2m a ω a √ . A simple application of the Heisenberg equation of motion leads to a set of equations of motion: and an adjoint equation: One can write these as one equation for a vector of four components using a Liouvillian operator L: dv/dt iLv, where v + ( a, b, a † , b † ) and Liouvillian: As is well known, one can diagonalize the Hamiltonian leading to an antisymmetric and a symmetric vibrational modes. Here this corresponds to diagonalizing the Liouvillian, Eq. 4, leading to the frequencies: where Ω 2 α 2 + (δω) 2 . We use as an example vibrational dynamics in the CS 2 molecule, where ω a ω b ω and Ω α. The eigenfrequencies, Eq. 5, become ω ± α/2 for the a ± b eigenvectors and similarly for their complex conjugate analogs. Taking expectation values it is useful to note that equations imply that the mean values at the time t can be computed given their initial values, for example: The mean values of the bond displacement, can thereby be computed. To describe the dispersion of the nuclear wave packets in coordinate space and in momentum space we need to extend our algebra. In the basis of the vibrational states the dispersion in the displacement of a local oscillator can be written as follows: And similar for the dispersion in the momentum space: Therefore we need to describe the dynamics of the { a 2 , ( a † ) 2 , a † a} † b} for the dispersion in R and P of the two coupled local bonds. Note connection between the dispersion of the wave packets in different representations to the energy transfer between the two coupled oscillators, 〈 a † a〉 or 〈 b † b〉 mean values. Therefore we need a hardware computing device that can emulate 10 logic variables. The computing hardware described above can emulate 16 logic variables when only the lowest four exciton bands are addressed and 64 when shorter laser pulses are used.

RESULTS
We use the device as discussed in section above to simulate the time dependence of the physical variables as revealed by the algebras for the coupled vibrations that are discussed in section III. Each coherence of the device is a point on the 2D frequency map generated by the 2D spectroscopy. The size dispersion of the quantum dots means that each point is actually a cloud of points at nearby frequencies. This dispersion in the frequency associated with each coherence is what enables one device to simulate different but similar physical systems. In this section we show how different coherences can simulate the time evolution of different physical variables of the system. The results for the mean bond displacement 〈R〉(t) and its dispersion σ 2 R (t), Eq. 7, are shown in panels (a) and (c) of Figure 3. The computations are for a symmetric system that is not stationary in time because of an asymmetry in the initial conditions. One oscillator starts with low energy and the other, the one shown, is initially energy richer. On the right side of the figure we show relevant segments of 2D frequency maps computed at different values of the time interval T that enable the device to simulate the two functions of time 〈R〉(t) and σ 2 R (t). The mean bond displacement of a local oscillator, frequency ω, in a symmetric molecule such as CO 2 , Figure 3A, is seen to be faster varying in time than the width, 3(c), of the wavepacket describing the oscillator. This reflects the difference in magnitude between the eigenfrequencies ω ± α/2 and the coupling parameter α. In the absence of coupling between the oscillators, α 0 in the Hamiltonian, Eq. 1, the width of a coherent state of the harmonic potential will not vary with time. The coupling induces a slow, frequency α, energy transfer between the two oscillators as reflected in the variation of the width of the wavepacket localized on one of the two bonds.
It might seem from Figure 3 that one has to assemble a new device for each different triatomic molecule. This is because for each one of the two frequencies ω + α/2 and α one would need a different electronic level structure in the dimer of the two quantum dots. The same size dispersion of the quantum dots that hitherto limited our abilities is here an advantage. The somewhat different sizes of the monomers means that there is a finite range of frequencies for each coherence. The range should not be too broad as otherwise we will lose resolution in frequency. Moreover, since the distribution of spacings of a particular coherence is, at low dispersions, Gaussian [3], we should not probe too far toward the wings. Complementary, when the range is broad we can only measure at short times before serious dephasing sets in. But the times cannot be too short because we want to measure only after the laser pulses are over. Practical implications of the distribution of spacings can be seen in Figure 4. Shown is a 2D frequency map for both regions used to compute the time-evolution of the bond distance (region A of Figure 4) and its dispersion (region B, Figure 4). Slightly different time-evolution of the coherences for the neighboring positions on the map enables fine-tuning of the period of the oscillations as is shown on the side panels of Figure 4.
A coherence is a complex valued observable so each coherence can describe two conjugate physical variables as shown in Figure 5. Equation 6 shows that the expectation values, 〈R〉(t) (〈 a〉(t) + 〈 a † 〉(t))/ 2mω √ and 〈P〉(t) i(〈 a † 〉(t)− 〈 a〉(t)) · mω √ / 2 √ of the bond displacement and of the momentum can both be computed at the same point of the map, by reading the real and imaginary values. The two values are shifted by a phase difference of π/2 as is to be expected. There is an uncertainty with the momentum, σ P , and position, σ R , as shown in Eqs. 7 and 8. Time-evolution of the uncertainty in the position, related to the width of the wavepacket in coordinate representation, is shown in Figures 3C,D. As discussed therein these uncertainties are time-dependent because of the coupling between the two oscillators so that they are more slowly varying and reflect the energy transfer between the two oscillators. The time-evolution of the uncertainty in the momentum representation is following the same time-dependence as σ R therefore no additional computation is needed.

CONCLUDING REMARKS
A versatile quantum mechanical computing device that operates on a laser addressed solid array of quantum dots dimers has been discussed. Experimental data shows that such a device can operate at room temperature [23,24]. The device was effectively used to simulate the quantal dynamics of nonstationary states of coupled vibrations. We show how the read-out of the coherences engineered in the device following interaction with a sequence of the laser pulses enables parallel computation of the two conjugated variables, mean position and momentum of the local oscillators. Extension to the more detailed dynamics is shown by simulating also the dispersion in coordinate representation of the quantum wave packets of the coupled oscillators. Simulation is provided by matching the beating frequency of the coherence and the oscillation frequency of the physical variable. This mapping is in principle always possible for systems where the algebraic description of the physical variables is closed upon their commutation with Hamiltonian. Diagonalization of the respective Liouvillian operator, Eq. 4, determines eigenfrequencies ω k , for example Eq. 5, and enables transformation to observables A k (t) that evolve in time as A k (0)exp(iω k t). It is these frequencies that need to be measured by the time evolution of the coherences of the device. The inevitable size dispersion of the quantum dots limits the span of time that is available for emulation but on the other hand allows a fine-tuning of a frequency of interest.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

AUTHOR CONTRIBUTIONS
RDL and FR designed the research; KG, RDL and FR performed the research and wrote the ms, KG, HG and FR carried out the computations.  Figure S4). The dispersion in sizes means that there is a range of frequencies corresponding to each coherence. Point A. A higher frequency coherence, more off the diagonal in the map, used to compute the time-evolution of the bond distance, Figures 3A,B. Shown are the intensities of the coherence vs time T for points in the range shown on the map. Point B. Coherence at a lower frequency for the computation of the bond dispersion, Figures 3C,D. A wider range in the frequency is possible but reading nearer to the diagonal is experimentally more challenging. Frontiers in Physics | www.frontiersin.org October 2020 | Volume 8 | Article 590699