Controlling the dynamics of ultracold polar molecules in optical tweezers

Ultracold molecules trapped in optical tweezers show great promise for the implementation of quantum technologies and precision measurements. We study a prototypical scenario where two interacting polar molecules placed in separate traps are controlled using an external electric field. This, for instance, enables a quantum computing scheme in which the rotational structure is used to encode the qubit states. We estimate the typical operation timescales needed for state engineering to be in the range of few microseconds. We further underline the important role of the spatial structure of the two-body states, with the potential for significant gate speedup employing trap-induced resonances.


I. INTRODUCTION
In recent years, ultracold polar molecules have been a subject of intense experimental and theoretical studies. Due to their rich internal structure and comparatively strong intermolecular interactions, they have been suggested as great candidates for quantum simulation [1][2][3] and computation [4][5][6][7][8][9], precision measurements of fundamental constants [10,11], as well as controlled chemistry [12][13][14]. The complex internal structure of molecules offers broad prospects for experimental control with external fields [15][16][17][18][19][20][21], but also leads to problems with cooling the system to reach quantum degeneracy. These include the lack of suitable cycling transitions (with several notable exceptions), and high inelastic collision rates leading to losses. Nowadays, after a series of experimental breakthroughs, an increasing number of groups can produce large ultracold molecular samples with high phase-space density, as well as trap them in an optical lattice or tweezer array with high filling and low entropy [22][23][24][25][26][27][28][29][30].
While general working principles of quantum engineering with trapped polar molecules are intuitive and based on well-established analogues in atomic systems [31][32][33][34], to provide meaningful experimental predictions it is essential to account for corrections with respect to commonly used approximations. One crucial aspect is the dipolar interaction, which is not only long-ranged but also state-dependent. As the rotational and hyperfine internal states strongly couple to electromagnetic fields, a suitable choice of molecular species and external conditions can lead to the realization of very diverse and rich many-body models [18,35]. This can be combined with strong optical confinement, which can also be made state-dependent due to the polarizability anisotropy of molecules which determines the trapping frequencies in an optical trap [36]. The characteristic interaction and confinement length scales compete with each other, making tight-binding and pseudopotential approximations questionable. As a result, even for a two-body problem, the full numerical solution requires extensive computational effort [37][38][39]. Precise information about the structure of states resulting from the interplay of strong interactions and confinement can nevertheless be very beneficial, as it allows to make use of the specific properties of the spectrum to increase the efficiency of state preparation and gate operation. One notable example are the so-called trap-induced resonances [40][41][42] resulting from the anticrossing between the molecular-like bound states and the spatially extended trap states.
In this work, we study the dynamics of a pair of ultracold polar molecules trapped in separate optical tweezers, fully taking into account the trap structure, internal rotational states, dipolar interaction, and external electric field. Our model is generic because it does not rely on any specific feature of particular molecular species and can be described in terms of a few characteristic length and energy scales. We study the system's evolution after a quench of the electric field value and set the stage for future calculations of the dynamics under optimized field pulses, providing estimates for the characteristic time scales.
The paper is structured as follows. In Section II we introduce the system Hamiltonian and discuss its properties relevant for state engineering. Then in Section III we analyze the spectra and dynamics and introduce the gate protocol, which is further discussed in Sec. IV. Conclusions are drawn in Sec. V. We provide the code supporting the findings of this study in Ref. [43].

A. System Hamiltonian
We consider a system of two polar molecules confined in separate three-dimensional potential wells. For simplicity, we assume the traps to be anisotropic harmonic, but anharmonic corrections are straightforward to include within our approach by extending the basis set to account for the coupled dynamics of the center of mass and relative motion as shown below. The molecules interact with each other via dipole-dipole forces, which are controlled with an external electric field.
Let us begin with a short discussion of the properties of a single polar molecule. Here we solely focus on its rotational structure within the rigid-rotor approximation and assume it remains in the ground electronic and vibrational state due to separation between the internal energy scales. In this work, we also neglect any possible hyperfine structure, as it is not relevant for the current demonstration. We choose the z axis along the electric field direction, leading to a simple description of the internal state where B is the rotational constant, j is the rotational angular momentum operator of a single molecule, E denotes the electric field magnitude, and d 0 =ê 0 · d is the z component of the electric dipole moment operator. In general, it is convenient here to use spherical tensor formalism where d p =ê p · d = d C 1 p (θ, φ), with d being the value of the permanent dipole moment in the body frame. The functions C k p denote the unnormalized spherical harmonics C k p (θ, φ) = 4π/(2k + 1)Y k p (θ, φ), whileê ±1 = ∓(x ± iŷ)/ √ 2 andê 0 =ẑ. When E = 0, the rotational angular momentum j is a good quantum number. Pure rotational states |jm have no mean dipole moment. At nonzero field only the projection of j onto the field axis, m, remains conserved and the eigenstates can be decomposed as jm = j c jj |j m , where jm denotes the state adiabatically connected to |jm . Still, j z jm = m jm . We will be mainly interested in the lowest lying rotational states connected to j = 0 and j = 1. Two natural parameters can be introduced here, the dimensionless field strength β = dE/B and the characteristic radius Let us now move to the case of two molecules. The total rotational angular momentum of the system is J = j 1 + j 2 with the projection onto the external field axis, M = m 1 + m 2 . The dipole-dipole interactions between two molecules are commonly represented as where rê r is the vector connecting the two molecules. This expression can also be conveniently rewritten by means of spherical tensors [44] which allows for separation of the part conserving the angular momentum projection The exchange term in this interaction potential can naturally flip an excitation or entangle the particles, being a starting point for a quantum gate scheme [6]. Here we instead rely on the state-dependent shift of the eigenstates in the presence of an external electric field for the purpose of quantum state engineering. We now discuss the spatial geometry of the system, which can be adjusted by changing the trap alignment with respect to the field direction and tuning the distance between the tweezer potential minima. Each optical tweezer is modeled as a cylindrically-symmetric anisotropic harmonic trap with frequencies ω z = ω, ω x = ω y = ω ⊥ = ηω (note that the z direction is chosen parallel to the electric field) with large anisotropy η. One can define here the characteristic trap length scales a ho = /µω z and l ⊥ = /µω ⊥ , where µ is the reduced mass of the system. For a harmonic trapping potential, the center-of-mass and relative motion are decoupled, and in our analysis we can focus on the relative motion described by the following Hamiltonian: where a is the separation between the two traps. For large anisotropy η 1 and sufficiently separated molecules, the transverse excitation is the highest energy scale in the system and the motion becomes effectively one-dimensional. In this regime, one can integrate out perpendicular degrees of freedom [45,46], assuming that the transverse wave function corresponds to the ground state of the harmonic trap. This yields an effective full one-dimensional Hamiltonian where with u = z/l ⊥ = zη 1/2 /a ho and erfc being the error function. Note that the part of the interaction involving internal states remains unchanged in the effective intermolecular interaction. Here, H rot denotes the internal state Hamiltonian (Eq. (1)) of the k-th molecule. In this work we followed the assumption that the transverse wave function is limited to the lowest oscillator mode, but expanding it on excited modes for a less anisotropic system is straightforward as we have checked that also for higher modes the coefficient in front of the delta function remains 8/3. If necessary, one can also add an additional short-range interaction to the effective potential in order to reproduce some physical scattering length.

B. Diagonalization
As a convenient basis for diagonalization of the Hamiltonian we take the states |i ≡ n z , j i where n z denotes the eigenstate of the harmonic oscillator in the z direction and j i 1(2) , m (i) 1(2) denote the rotational states of the first (second) molecule. The function |n z can be centered either at z = 0 or at z = a. The former choice is better suited to describe the case in which the molecules lie close to each other or even form a bound state, while the latter should work well if the interaction is weak and the particles are well separated. We checked that at the length scales relevant for our case, better numerical stability is achieved using states centered at z = 0, while at larger separation or weaker interactions, it would be beneficial to use the other basis.
Let us now briefly discuss the matrix elements i| H |i of the Hamiltonian of Eq. (6). Starting with the trapping potential along the z-axis, H trap of Eq. (5), the integral is calculated analytically, while the spatial part of the interaction potential, i| V eff dd |i , is calculated numerically. The matrix element of the rotational angular momentum for the k-th molecule are  The effect of the electric field on the k-th molecule and the matrix elements of the dipolar interaction can be calculated from the definition of the dipole moment operator [44]. For completeness, we provide here the matrix element of d q required to calculate both quantities where (a, m a ; b, m b |c, m c ) are the Clebsch-Gordan coefficients.

A. Energy spectrum
We now discuss the properties of the energy spectrum of the system. For the physically realistic case in which ω B, the eigenstates separate into branches corresponding to different numbers of rotational excitations. As the most intuitive experimental control knobs for the system are the trap separation and electric field magnitude, we study the energy levels as a function of these parameters. For our numerical calculations we use the values corresponding to NaCs molecules with d = 4.607 Debye and B = 1.813 GHz, taking the trap frequency ω = 2π · 50 kHz and η = 10. This implies a ho ≈ 960.7 a 0 .
First, we investigate the role of trap separation in figure 2 which shows the three lowest branches corresponding to the total angular momentum projections M = m 1 + m 2 = 0, 1, and 2. In all cases, one can clearly distinguish two types of eigenstates: the trap states, whose energy is roughly independent of the distance a; and the bound states for which the energy goes up roughly as a 2 . This behavior is typical for the chosen system geometry [47] and does not depend on the specific type of interactions. Here, in contrast to the commonly studied contact interaction case, strong attractive potential well leads to the emergence of multiple bound states and large energy shifts. Suitable states for quantum gate realization correspond to the trap states that can be efficiently prepared in remote traps and then brought together. Note that the dipolar bound states display anticrossings with the trap states. This phenomenon, called the trap-induced resonance [40], stems from interaction-induced coupling and in the energy spectrum looks similar to a Feshbach resonance. It can be utilized for different purposes such as production of molecular states, but also to shape the energy spectrum by shifting the energy of a trap state in one of the branches while leaving the other intact. For example, in figure 2b there are more than five resonances with the lowest trap state, while only two are visible in panel 2a. At small distances, all types of states become strongly mixed.
It is important to note that the second branch corresponding to the total angular momentum projection M = 1 experiences much stronger effects of the dipole-dipole interaction than other branches. It can be understood on the basis of perturbation theory calculations, where the dipole-dipole interaction has a non-zero effect in the first order only for pairs of states with j 1 = 0(1), j 2 = 1(0), while for j 1 = j 2 = 0 and j 1 = j 2 = 1 it contributes only as the second-order correction.
Turning on the electric field induces a nonzero net dipole moment in the molecules and thereby leads to stronger interactions. This is reflected in the spectra as strong shifts of the energy of all two-body states to lower values, as illustrated in figure 3, which shows the energies of the lowest lying trap states in different branches. As we have chosen the field direction to be parallel to the trap axis, the interactions are attractive. Here we chose the trapping potential to be the weakest along the dipole axis as well as the trap alignment, which enhances the attraction and simplifies the spectrum, allowing us to neglect the transverse excitations. In general, it is more convenient experimentally to realize a different setting, while our setup would require an additional light sheet or a lattice. For the general case the density of levels would increase and the spectra would become more complicated, especially the trap-excited states. However, the trap-induced resonances would still be present.

B. Dynamics
Having understood the basic features of the energy spectrum, we now proceed with the dynamics. We note that recent proposals for quantum engineering protocols involving polar molecules include taking advantage of the dipolar exchange interaction [6] or utilizing a microwave pulse [7]. Our scheme is complementary to these approaches, being based on applying an electric field pulse. The interaction with the field leads to state-dependent energy shifts, allowing in principle for the realization of various quantum gates, but optimization of the pulse shape would be necessary to achieve the desired phase accumulation. Here we will focus our analysis mainly on the simple scenario of electric field quenches and study the population of motionally excited states during the process. While such excitations can be regarded as a fidelity leakage source, it is possible to take them into account and design suitable control pulses that will not only keep the final state close to the trap ground state but also use the full space of states as a resource for gate speedup [48,49].
In order to perform a quantum gate, one has first to choose the suitable computational basis. States with varying M are natural candidates for this purpose, as they are not coupled with each other by the dominating interaction term of equation (4) that we consider here. We thus choose to focus on the two-particle states that have the largest overlap with the pair of molecules being trapped in the motional ground state of separated tweezers. The single qubit states can then be chosen as |0 = |j = 0, m = 0 and |1 = |j = 1, m = 1 . Rotational excitations are long-lived and thus very suitable for our purpose. Then the two-qubit basis is composed of states |00 with J = 0, M = 0, |+ being the symmetric combination with J = 1, M = 1, |− being the antisymmetric state from the same branch and |11 with J = 2, M = 2. In each case we assume the state is initially prepared in a motional eigenstate with the corresponding amplitude c i (i for initially occupied). The evolution can lead to occupation of multiple other trap states which we will denote with the c f coefficient (f for final).
We start the calculation by diagonalizing the Hamiltonian for a zero electric field, as in our scenario of interest the field strength is the only parameter that varies during the process. This solution has the advantage of being conceptually simple, while providing short operation times. Local manipulation of the qubits can be achieved by individually addressing the molecules e.g. with an off-resonant laser. Then we move to the interaction picture and solve the corresponding Schrödinger equation numerically (using the solver implemented in Mathematica) for the given time dependence of β(t). As before, we choose NaCs molecules as an example and set the distance between the traps to 10.4 a ho , which for the system parameters we consider equals 10 4 a 0 , away from the trap-induced resonances in the M = 0 and M = 2 branches visible in figure 2.

Time evolution with constant electric field
We first study the quench scenario in which the electric field is suddenly turned on and remains constant throughout the evolution. The operation time is set to 200 ns. We will see that this time is long enough for multiple oscillations of the wave function components.
Let us start with the case of a weak field, shown in figure 4 for β = 2.3 · 10 −2 . The evolution of M = 0 and M = 2 states in panels (a) and (c) shows similar oscillatory behavior, while the M = 1 state in panel (b) undergoes more complicated dynamics. This can be explained by the larger density of states in this branch due to the strong dipolar attraction which mixes trap and bound states. In general, for the electric field values studied here, we observe that the evolution does not lead to coupling of the initially occupied trap state to the bound states and does not excite higher rotational branches. The main couplings occur between the nearby trap states which have the largest overlap with the starting one. For the case depicted in figure 5 the field is increased to β = 0.16. Here the evolution becomes much faster and the population spreads over a larger number of motional basis states. This is once more especially visible for the M = 1 state, where the population of the initially occupied state drops to below 20% as shown in panel 5(d). For reference we also show in panels 5(e), (f) the occupation of two states that are initially not present, while their population arises due to strong overlap with the initial state after the quench.

Shaping the pulse
While the quench scenario is instructive, the ultimate goal would be to switch the field on and off in a smooth way and control its shape in order to achieve desired operations. In order to realize a quantum gate, it is required that the internal levels acquire specific state-dependent phases (e.g. a CNOT operation), while the final motional state should not differ from the initial one. As an initial step towards the full quantum computing proposal, here we study the evolution of the system under a simple pulse β(t) = β 0 C(t) sin πt τ with τ being the operation time, β 0 being the pulse amplitude, and C(t) being a correction written as truncated Fourier series C(t) ∝ i A i cos(ξ i t) + B i sin(ξ i t) in the spirit of the chopped random basis optimization method [49,50]. We restrict the total operation time to 150 ns for this simple demonstration. For an estimate of the quantum speed limit (the time required for achieving close to unit fidelity), one can look at the inverse of the smallest energy level separation. This is roughly given by the trap frequency, which here corresponds to 50kHz leading to ∼ 20 µs gate times.
The target of the operation that we chose for the demonstration was to perform a controlled phase gate (in the present case realized in such a way that three states acquire a π phase, while the |11 state does not) in the qubit space with only a few optimization steps using gradient search in a small basis set. The system evolution for this case is shown in figure 6. One can notice that if the trap excitations were disregarded, the pulse would already reach 84% fidelity after performing only two optimization steps. However, the |11 state corresponding to M = 2 becomes transferred to an excited motional state as a result of the evolution (see figure 6(d) and (e)), such that the actual fidelity is zero. This shows that more elaborate control schemes would need to be applied in a realistic system as the motional decoherence can be a problem (note that neglecting the trap dynamics or assuming a gaussian spatial profile with some finite thermal width disregards this problem). The evolution of the |00 state in panel 6(a) is notably slower than the others, as its energy is closest to zero, resulting in low oscillation frequency in the interaction picture. For the higher rotational branches the states chosen as the computational basis are more excited and undergo a complicated  evolution. The pulse shape shown in panel 6(f) is slowly varying and thus experimentally realistic. More extensive calculations using larger basis sets, more optimization steps, and realistic trapping potentials will follow in future work.

IV. DISCUSSION
The numerical examples above indicate that while the dynamics does not lead to couplings between the states with different rotational quantum numbers, the trap states can become mixed during the operation. This would lead to the spreading of the wavepacket and limit the number of possible operations. On the other hand, quantum control techniques can be applied to ensure high fidelity [48]. The main question for the system at hand is not the possibility of realizing quantum gates but rather the quantum speed limit achievable in various settings [51]. In a generic setup of two trapped molecules, neglecting possible experimental imperfections, the most important source of fidelity loss lies in the strong dipolar couplings between the motional states from the same rotational branch. This, in principle, affects any possible gate scenario that relies on the dipolar interaction.
Interestingly, the trap-induced resonance phenomenon shown in figure 2 can be utilized to manipulate the energy spectrum and increase the gap between the chosen qubit state and other trap states. Close to the resonance, the trap acquires some bound state character such that part of the wave function is localized at small interparticle separation, while the other states remain delocalized over the trap wqhich reduces their coupling. Furthermore, the energy shift resulting from the anticrossing will be translated to an additional phase shift of the affected state with respect to another qubit state, providing the possibility for speedup.

V. CONCLUSIONS
We have analyzed the prospects for quantum state engineering of ultracold polar molecules trapped in separate optical tweezers and controlled using an external electric field. By taking into account the complete structure of the trap states, we have shown how motional excitations could arise during the evolution. This allows for the implementation of more elaborate control schemes, which in principle would allow for working close to the quantum  speed limit set by the harmonic confinement frequency. Reaching this goal can be made easier by utilizing trap-induced resonances, which strongly depend on the internal state of the molecules and thus can be precisely tuned. Possible extensions of the present work include a study of the resonances in experimentally realistic traps and taking into account more details of the interaction potential to deliver precise predictions on the trap-induced resonance positions. Then, optimal control techniques can be utilized to design fast gate protocols with high fidelity. In addition, we suppose that considering more details of particular molecular species, such as including their hyperfine structure and adding more external fields (e.g., microwave) to optimize the qubit space further, will ultimately lead to a full quantum computation toolbox, as well as allow for more detailed studies of molecular interactions.