Phonon arithmetic in a trapped ion system

Single-quantum level operations are important tools to manipulate a quantum state. Annihilation or creation of single particles translates a quantum state to another by adding or subtracting a particle, depending on how many are already in the given state. The operations are probabilistic and the success rate has yet been low in their experimental realization. Here we experimentally demonstrate (near) deterministic addition and subtraction of a bosonic particle, in particular a phonon of ionic motion in a harmonic potential. We realize the operations by coupling phonons to an auxiliary two-level system and applying transitionless adiabatic passage. We show handy repetition of the operations on various initial states and demonstrate by the reconstruction of the density matrices that the operations preserve coherences. We observe the transformation of a classical state to a highly non-classical one and a Gaussian state to a non-Gaussian one by applying a sequence of operations deterministically.

I n quantum mechanics, bosonic creationâ y and annihilationâ operators bear the following operator relationŝ where n j i stands for a Fock state of n bosons. The proportionality factors ffiffiffiffiffiffiffiffiffiffi ffi n þ 1 p and ffiffiffi n p appear due to the symmetric indistinguishable nature of bosons 1,2 . Thus, the addition or subtraction in quantum domain involves the modification of the probability amplitude of state due to the excitation n-dependent factor. The bosonic annihilation and creation operation have been proposed to be a building block to generate an arbitrary quantum state 3,4 , to distill entanglement or non-locality 5 and to transform a classical state to a nonclassical state 6 . In recent times, there have been seminal works to realize the bosonic operations at the single-boson level for the test of foundations and applications of quantum mechanics [7][8][9][10][11][12][13] . However, such bosonic operations are not trace preserving and inherently probabilistic. The probability of success has so far been extremely low in their implementation for photonic fields. A higher probability may be obtained at the expense of lowering the performance fidelity 14 . Owing to this, a repetition of these operations has been limited.
The conventional addition and subtraction of a particle can be written asŜ The additionŜ þ takes the n particle state to the (n þ 1) state, whereas the subtraction operationŜ À brings the n state to the (n À 1) state without incurring additional factors. Seeing the form of the operations in equation (2), we immediately recognize that S þ is a deterministic process whileŜ À may not be, as it is not possible to subtract a particle from vacuum |0i. When the vacuum component of the initial state is small, the subtraction can be done near-deterministically. In recent times, there have been theoretical proposals of the operations (equation (2)) for the generation of an arbitrary quantum state 15 , the measurement of vacuum 16 , the transformation to a non-classical state 17 and the amplification of a quantum state 18 . In particular, such arithmetic operations form an important component of a qubit gate operation for ions in a harmonic potential 19 . The operations (equation (2)) were also suggested as the elements of a phase operator 20 . Beyond the quantum-state engineering, the subtraction can be used for the sub-Doppler cooling in the trapped ion system 21 .
In this study, we experimentally demonstrate deterministic addition and near-deterministic subtraction of a bosonic particle, in particular a phonon of a 171 Yb þ ion trapped in a harmonic potential. We realize the operations by coupling phonons to an auxiliary two-level system, so called the hybrid scheme of discrete and continuous variable 22 , and applying a transitionless quantum-driving scheme. We perform the operations on superpositions of Fock states and coherent states, and we demonstrate that our single-phonon operations are (near) deterministic and preserve coherence. By applying a sequence of operations deterministically, we show that classical states are turned into non-classical ones manifesting highly sub-Poissonian statistics and negativity in the Wigner function.

Results
Experimental implementations. We implement theŜ þ andŜ À operations of equation (2) on a vibrational mode of frequency o X for a single trapped 171 Yb þ ion in a three-dimensional harmonic potential 23 through its interaction with the two-level system of atomic energy levels. The harmonic potential is generated by an oscillating electric field in the radial axis with trap frequency o X ¼ (2p)2.8 MHz. The two-level system is represented by two hyperfine states F¼1; m F ¼0 j i " j i and F¼0; m F ¼0 j i # j i of the S 1/2 manifold with the transition frequency o HF ¼ (2p)12.6428 GHz. As shown in Fig. 1a, the anti-Jaynes-Cummings (aJC) interaction or blue-sideband transition, wŝ þ e iDt þ h:c:, is realized by the stimulated Raman laser beams with beat-note frequency Here, O is the Rabi frequency of the two-level system, Z¼ Dk ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi ' =2Mo X p the Lamb-Dicke parameter, Dk the net wave vector of the Raman laser beams and M the mass of 171 Yb þ ion. The aJC coupling produces the transition between #; n j i and "; n þ 1 j iwith the oscillation frequency of ffiffiffiffiffiffiffiffiffiffi ffi n þ 1 p ZO, where the ffiffiffiffiffiffiffiffiffiffi ffi n þ 1 p factor comes from the fundamental property ofâ w andâ operators in equation (1). Therefore, the application of the simple aJC interaction does not transfer #; n j i to "; n þ 1 j i in an n-independent manner at a fixed duration of time.
The full population transfer independent of the initial motion state, the uniform blue-sideband transition, P n¼0 "; n þ 1 j i #; n h jþ h:c:, can be obtained by the application of the stimulated Raman adiabatic passage [24][25][26] . The adiabatic scheme provides a robust transfer against the variations in the transition rate either from the intensity change of the control Raman beams or from the property of the transitions 25 . Therefore, a properly designed adiabatic passage would allow a decent state transfer for a wide range of phonon number states through the aJC interaction, despite the ffiffiffiffiffiffiffiffiffiffi ffi n þ 1 p ZO dependence. In the adiabatic passage, typically ZO slowly increases at the beginning and decreases at the end, that is, O(t) ¼ O 0 sin(pt/T), whereas the detuning D changes according to D(t) ¼ D 0 cos(pt/T), where T is the total transfer time, across the resonance. However, for the applicability of the scheme to a wide range of initial phonon numbers with high fidelity, we should set D 0 as high as ffiffiffiffiffiffiffiffiffiffiffiffiffi ffi n M þ 1 p ZO 0 =2, where n M is the largest phonon number for the transfer and should fulfill the adiabatic condition, T ) 1=ZO 0 . In our experimental conditions, the reasonable duration T for such an adiabatic transfer is around 21 times that of p-pulse duration for the blue-sideband transition of the ground state p/ZO 0 , when we include up to the maximum phonon number n M ¼ 6.
In recent times, the transitionless quantum-driving scheme [27][28][29][30] has been developed to speed up the adiabatic control. When the transfer is non-adiabatic, it introduces an additional term in the Hamiltonian of the instantaneous basis, which is wŝ þ e iDt þ h:c:, where b is in the order of p/ZO 0 /2T, the ratio between the p-pulse duration and the total operation time. By adding a counter-diabatic term in the control, we can suppress the non-adiabatic excitation with a reasonable speed up over the adiabatic passage. For the aJC interaction, the optimal values of D 0 and b are dependent on the phonon number n for the given O 0 . In our experiment, we optimize D 0 and b for the case of the geometric average of the minimum and maximum phonon number, . By doing this, we are able to reduce the total duration of the operation from 21 to 7 times p/ZO 0 without sacrificing the fidelity of the rapid adiabatic passage for the same range of phonons n M ¼ 6.
To rigorously achieve the uniform blue-sideband operation P n¼0 "; n þ 1 j i#; n h jþ h:c:, it is important not only to increase the phonon number but also to preserve the relative phases between component states of the quantum state. For our previous realization 23 , the different extra phases were accumulated depending on the phonon number of the initial state, which prevented from keeping the initial phase coherences (see 'AC Stark shift compensation' in Methods). In this study, we have ...  ARTICLE developed a sequence of operations compensating the phononnumber-dependent phases based on the spin-echo principle. As shown in Fig. 1b, we invert the sign of O and reverse the control of D in the middle of the sequence of the operation, which symmetrizes the whole operation and produces the accumulated phases of opposite signs before and after the inversion and reverse. Therefore, the total phases are cancelled out at the end of the operation (see 'AC Stark shift compensation' in Methods).
Conventional phase-coherent addition. We implement the coherent addition operationŜ þ in equation (2) by applying the uniform blue-sideband transfer followed by the p-pulse of carrier transition as shown in Fig. 2a. By this operation, we add single phonons deterministically, independent of the initial number of phonons. To manifest the coherence property of the addition operations, we prepare an initial state 1 ffiffi 2 p n¼0 j iþ n¼1 j i ð Þ , apply the operations up to three times and measure the density matrices of the resulting phonon states. As shown in Fig. 2b, the coherences represented by the off-diagonal terms of the density matrix clearly remain after the multiple addition processes up to three times.
We also consider an initial coherent state of the amplitude a ¼ 0.81, which is achieved by displacing the vacuum. The density matrix is reconstructed (see Fig. 2c,d)   distributions for eight different angles (see 'Reconstruction of the phonon density matrix' in Methods). The phonon number distributions are obtained by observing the time evolutions of the standard blue-sideband transitions, similar to the direct reconstruction scheme of the phonon density matrix 32 . One immediate consequence of adding single phonons on the coherent state is the production of sub-Poissonian phonon statistics, because the addition increases the average phonon number but not the shape of the distribution nor the variance. We observe that the ratios between the variance and the average phonon number reduce from 1 to 0.43, 0.39 and 0.2, after 1-3 single-phonon additions, respectively. Applying the first addition operation, we detect negativity in the Wigner function as shown in the second column of Fig. 2d. It is important to note that the addition operation, which converts a coherent state to a highly non-Gaussian and non-classical state, is deterministic unlike the case ofâ y operation 9-11 . There is a technical limit in the number of additions we can apply, owing to the validity of the adiabatic approximation and the heating process of phonons 33 . Under this limitation, we could perform the operations three times without the significant loss of fidelity. As shown in Fig. 2d, the experimental results and the theoretical predictions for the Wigner functions are in agreement, except the shift of the orientation, which would come from the phase offset between the preparation and the analysing pulses. This is significant in comparison with the photonic realization of bosonic operations of single-photon creation and annihilation 11 . We detect negativity in the Wigner function as shown in the second column of Fig. 2, where we obtain the Wigner function of the state from the reconstructed density matrix.
Conventional subtraction operation. The subtraction operation S À in equation (2) is realized by reversing the sequence of the addition operation, that is, the application of the p-pulse of carrier transition and the uniform blue-sideband transfer, followed by the fluorescence detection as shown in Fig. 3a. The zero phonon state n¼0 j i is eliminated after the subtraction, which is implemented by the conditional measurement in our experimental scheme. After the detection sequence, we only collect the data with no fluorescence, which has the success rate given by the probability of the non-zero phonon states. We examine the performance of the subtraction operation with an initial phonon superposition state 1 ffiffi 2 p n¼2 j iþ n¼3 j i ð Þ . As shown in Fig. 3b, the subtraction operation reduces the phonon excitation by one quanta. After the second application of the subtraction, the off-diagonal terms of the density matrix are significantly reduced, which shows the current limit in experiments due to heating of the system. We also prepare a coherent state a¼1:2 j i and apply the subtraction twice as shown in Fig. 3c, which shows qualitatively that the subtraction operation works for any initial quantum state. The subtraction operation can squeeze a coherent state, which is different from annihilationâ that has the coherence state as its eigenstate. However, our experimental precision is not high enough to observe the squeezing effect.
Commutation relation of addition and subtraction operations. We study experimentally how the quantum states are changed depending on the order of the addition and subtraction for an initial coherent state c i j i¼ a¼1:2 j i. If we add and then subtract S ÀŜ þ c i j i, the state after the sequence is the same as the original one, as there is no amplitude modification. For the case of subtraction-then-additionŜ þŜ À c i j i the final state does not have vacuum component, because a vacuum component is annihilated byŜ À and cannot be restored byŜ þ . Figure 4b shows the experimental result ofŜ ÀŜ þ c i j i, which is basically identical to the initial state of Fig. 4a. Figure 4c shows the result after the operation ofŜ þŜ À c i j i, where there is no significant vacuum component in the density matrix. The vacuum component is not perfectly absent because of the detection error during the projective measurement based on the atomic fluorescence and heating of the system. The fluorescence detection duration is comparable to the motional coherence time of our system, which makes the off-diagonal part of the final state suppressed significantly (see 'Error analysis' in Methods). Our experimental result is well in line with the non-commuting relation of the Susskind-Glogower's phase operators, that is,

Discussion
Owing to the capability of deterministically generating a nonclassical and non-Gaussian states, the conventional addition and subtraction operations provide an efficient scheme for quantum state engineering 34,35 . The current work can also be a stepping stone to realize the ion qubit gate operation as proposed by Schneider et al. 19 without having to worry about long operation times and loss of coherences. This scheme may be further applied to various quantum optics setups such as cavity quantum electrodynamics and optomechanics. We note a theoretical scheme that the non-Gaussian information of a phonon state would be transferred to that of photon through phonon-photon coupling 36,37 . It has also been discussed that the conventional arithmetic addition and subtraction can be used to measure the vacuum state without disturbing the state 16 and therefore it can ARTICLE be directly applied to construct the Q-function of a quantum state of motion. We have noted that an arithmetic subtraction was performed for a photonic system 38 .

Methods
AC Stark shift compensation. The AC Stark shift in the adiabatic operations mainly comes from the off-resonant coupling to the carrier transition, the transition between S 1/2 2P 1/2 states of 171 Yb þ ion and the other radial motional mode (o Y Eo X þ (2p)0.4 MHz). The dominant AC stark shift comes from the carrier transition of frequency 38.5 kHz and the Lamb-Dicke parameter Z ¼ 0.089. The amount of the shift brought by the Y mode is given by O 2 0Y 2 oY À oX ð Þ , that is, B20 times smaller than that from the carrier coupling. The AC stark shift between qubit states from the Raman laser beams due to S 1/2 2P 1/2 transition is oHF DR $ 2p ð Þ1 kHz, where g R1 and g R2 are the coupling strengths of Raman 1 and Raman 2 beams, respectively, and the D R ¼ (2p) 18 THz is the detuning from the level of 2 P 1/2 .
We consider total AC stark shift as the form of O t ð Þ j j 2 2Dtotal , where D total is the detuning effectively including all the possible origin of AC stark shifts discussed above. We obtain o . The actual frequency of blue-sideband o act bsb is measured by observing the resonant excitation. Including the AC stark shift, the actual waveform of O(t) that we apply on our arbitrary waveform generator is as follows: Here, D(t 0 ) ¼ D 0 cos(pt/T) and we note that the imaginary part in original form of is changed to sin-wave, which has the p 2 phase difference. Here, f(t) is calculated as follows: Reconstruction of the phonon density matrix. We use an iterative algorithm in ref. 31 for the reconstruction of a quantum state. It consists of a maximumlikelihood estimation solved by the expectation-maximization algorithm followed by a unitary transformation of the eigenbasis of the density matrix r. Based on the measured phonon distribution f n by N measurements after a displacement of the quantum state, we aim to get real probabilities p n ¼ n h jr n j i that are as close to the observed frequencies f n as possible, subject to the maximumlikelihood functional from which we reconstruct r. This likelihood functional can be interpreted as a linear and positive problem in the classical signal processing: where r i is the eigenvalues of r and h in is a positive kernel. We can solve this linear and positive problem with the expectation-maximization algorithm 39,40 : which is initially set to a positive vector r r i 40 8i ð Þ . This is repeated for different displacements.
The second part of the reconstruction scheme aims at getting the eigenbasis diagonalizing the density matrix. This part consists of two steps: reconstruction of the eigenvectors of r in a fixed basis and rotation of the basis using a unitary transformation with the infinitesimal form U e iEG % 1 þ iEG and E is a small positive real number. G ¼ i[r, R] is chosen as a Hermitian generator of the unitary transformation, where R is a semipositive definite Hermitian operator R¼ P n fn pn n j i n h j.
Starting from some positive initial density matrix r, we continue repetition of first finding eigenvalues r i using the expectation-maximization iterative algorithm (equation (9)) and then finding eigenvectors f i by unitarily transforming the old ones. The likelihood of the estimate p n is increased and we finally reach to determine the density matrix r ¼ X n r n f n j i f n h j: ð11Þ Error analysis. Dominant error comes from the phonon heating process caused by the electric-field noise from the trapped electrodes 41 . Heating decreases the Fock state preparation fidelity and affects the adiabatic blue sideband process. Its time evolution is known to be described by 41,42 : y r t ð Þâ À r t ð Þââ y Àââ y r t ð Þ þ g 2 n þ 1 ð Þ 2âr t ð Þâ y À r t ð Þââ y Àâ yâ r t ð Þ : where g is the coupling strength between the ion motion and the thermal reservoir and n is the average phonon number for the thermal reservoir. In our experimental setup, the heating rate g n % gn Ç 1 is 150 Hz. It can be reduced by using a large trap, cleaning the electrodes 43 (equivalent to reducing g) or cooling the trap (equivalent to reducing n) (ref. 44). Error may be caused by fitting, as some noise in the blue-sideband curves may be recognized as a high phonon population. Although the Fock state preparation error will be involved in the transition probability, fortunately the error for a small n is relatively insignificant and the population mainly resides in small n-values for the initial states we prepared in the experiments.
There are some other small errors that we did not calculate or simulate, including the fluctuation of the trap frequency B0.1 kHz coming from B2% intensity fluctuation of Raman laser and B1 kHz trap frequency fluctuation coming from sudden changes of the ion position. We see how much our experimental results are modified when comparing with the simulations 23 .